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Abstract 



The recent years have seen a tangible progress in the fermion field discretiza- 
tion with an exact chiral symmetry, which is realized on the lattice through 
the Ginsparg- Wilson relation. At present the simulation costs of Ginsparg- 
Wilson fermions are considerably higher than for the conventional formula- 
tions without the exact chiral symmetry and therefore only quenched sim- 
ulations in Quantum Chromodynamics (QCD) are feasible so far. However, 
the chiral fermions are especially important to address the ab initio QCD 
calculations at small quark masses of the low energy constants of the chiral 
Lagrangian. This thesis is dedicated to explore the feasibility of such calcu- 
lations in the e-regime of QCD for the extraction of physical information. 
We apply two formulations of the Ginsparg- Wilson fermions, namely, the 
Neuberger operator and the hypercube overlap operator to compute the ob- 
servables of interest. Since the calculations are done in the e-regime, each 
observable has to be sampled within a fixed topological sector, which is well 
defined for the Ginsparg- Wilson fermions as the index of the Dirac operator. 
As a main result we present the comparison of the distributions of the leading 
individual eigenvalues of the Neuberger operator in QCD and the analytical 
predictions of chiral random matrix theory. We observe a good agreement as 
long as each side of the physical volume exceeds about L pa 1.12 fm. At the 
same time the chiral condensate £ can also be estimated. It turns out that 
this bound for L is generic and sets the size of the physical volume where 
the axial correlator behaves according to chiral perturbation theory. This al- 
lows us to compute a value for the pion decay constant F w . The simulations 
also show that due to the high probability of the near-zero modes it is pro- 
hibitively difficult to sample the axial correlator in the neutral topological 
sector. In the higher sectors, however, we observe that the sensitivity of the 
analytical predictions for the axial correlator to extract E is lost to a large 
extent. As an alternative procedure we only consider the contribution from 
the zero modes. Here we are able to obtain an estimate for F n and a, where a 
is a low energy constant peculiar to quenching. We calculate the topological 
susceptibility, both for the Neuberger operator and for the overlap hypercube 
operator. It turns out that the result with the overlap hypercube operator 
is closer to the continuum limit. Also the locality properties are superior to 
those of the Neuberger fermions. As a theoretical development the Liischer 
topology conserving gauge action is investigated. This enables us to sample 
the observables of interest in the e-regime without recomputing the index. 
We can report that a promising gauge action has been identified. 



Zusammenfassung 



In den letzten Jahren gab es konkrete Fortschritte in der Diskretisierung 
von Fermionfeldern mit exakter chiraler Symmetrie. Diese Symmetrie wird 
auf dem Gitter durch die Ginsparg- Wilson Beziehung realisiert. Momen- 
tan sind alle bekannten Formulierungen wesentlich aufwandiger als kon- 
ventionelle Fermionen ohne chirale Symmetrie. Deswegen sind heutzutage 
nur ,, quenched" Simulationen der Quantenchromodynamik (QCD) moglich. 
Doch die chiralen Fermionen sind besonders wichtig fur ab initio Berechnun- 
gen der Niederenergie-Konstanten der chiralen Lagrangedichte bei kleinen 
Quarkmassen. Das Ziel dieser Dissertation besteht darin, die Realisierbarkeit 
dieser Methode zur Gewinnung physikalischer Informationen im e-Regime 
der QCD zu erforschen. Wir haben zwei Formulierungen eingesetzt. Dies sind 
der Neuberger Operator und Overlap Hyperkubus Operator. Da die Rech- 
nungen im e-Regime durchgefurt werden, muss man jede Observable in einem 
bestimmten topologischen Sektoren messen. Bei Ginsparg- Wilson Fermionen 
ist dieser durch den Index des Dirac Operators defmiert. Ein Hauptergebniss 
dieser Arbeit ist der Vergleich der Wahrscheinlichkeitsverteilungen einzel- 
ner Eigenwerte des Neuberger Operators in der QCD mit entsprechenden 
analytischen Vorhersagen der Theorie der Zufallsmatritzen. Wir beobacht- 
en eine gute Ubereinstimmung solange jede Seite L des physikalischen Vol- 
umens grofier als etwa 1.12 fm ist. Dabei kann auch das chirale Konden- 
sat X abgeschatzt werden. Es ergab sich, dass diese untere Schranke von 
L allgemein gilt und die GroBe des physikalischen Volumens, auf dem der 
Axialkorrelator den Vorhersagen des chiralen Storungstheorie folgt, festlegt. 
Damit konnen wir die Pionzerfallskonstante F n bestimmen. Unsere Simula- 
tionen zeigen, dass wegen der grofien Wahrscheinlichkeit niedriger Eigenwerte 
die Messung des Axialkorrelators im topologischen neutralen Sektor extrem 
aufwandig ist. Doch reicht die Empfindlichkeit der Vorhersagen der chiralen 
Storungstheorie in hoheren topologischen Sektoren bei der gegebenen Statis- 
tik nicht zur Bestimmung von £ aus. Als alternative Methode, gehen wir 
dazu iiber, allein den Beitrag der Nullmoden zu betrachten. Hier konnen 
wir Abschatzungen fur F n und a gewinnen. Dabei ist a eine Niederenergie- 
Konstante die nur in der ,, quenched" Naherung erscheint. Wir berechnen 
die topologische Suszeptibilitat fur den Neuberger und Overlap Hyperkubus 
Operator. Im letzten Fall ist der berechnete Wert naher beim Kontinuum- 
slimes. Die Lokalisierung fur den Overlap Hyperkubus Operator ist auch 
besser als fur den Neuberger Operator. Unser anderes Ziel ist die Erforschung 
einer topologieerhaltenden Eichwirkung. Hier konnen wir berichten, dass eine 
vielversprechende Eichwirkung gefunden worden ist. 
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Chapter 1 
Motivation 



This dissertation is dedicated to the Quantum Chromodynamics (QCD), a 
theory that explained already many effects and processes in the hadron world. 

QCD is a part of the Standard Model of elementary particles with the 
gauge group £77(3) © SU (2) © U(l). It is the theory that describes the strong 
interactions of quarks and gluons. Quarks are the fermion fields which belong 
to the fundamental representation of the gauge group £77(3) and gluons 
being in its Lie algebra. In QCD the strong interactions are mediated by the 
exchange of these gauge bosons. Due to the non-Abelian nature of the gauge 
group £77(3) a self interaction of the gauge bosons is also present. There is 
experimental evidence that the quarks come in six flavors: up (u), down (d), 
strange (s), charm (c), bottom (b) and top (t). 

Many QCD calculations were successful in the domain of high energies 
where — due to the asymptotic freedom — the coupling constant becomes a 
convenient small expansion parameter. This is the perturbative regime where 
many quantities of interest can be expressed in terms of the asymptotic 
expansion in powers of the strong coupling constant a s . 

In the domain of energies of order 1 GeV the coupling constant rises 
to order 1, 0(1), and the perturbative treatment fails. Therefore a non- 
perturbative approach must be applied to address the problems of interest in 
this domain. This scale of energies corresponds to color singlet bound states 
of quarks when they form hadrons. One distinguishes between mesons, con- 
sisting of a quark and antiquark, and baryons, consisting of three quarks. An 
elegant classification of the mesons and baryons is given by the quark model. 
The particles are identified with the fields belonging to the irreducible rep- 
resentation of their flavor group. To this end the classification of the mesons 
for the flavor group SU(2) is 2 © 2 = 3 © 1, i.e. there is a triplet of meson 
fields and a singlet. For the SU (3) group the pattern is 3 © 3 = 8 © 1, i.e. an 
octet and a singlet. In experiments on the deep inelastic scattering the build- 
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ing blocks of hadrons, the partons, were identified. The quarks are, however, 
observed in nature only in such bound states and no free quark has ever been 
found. This phenomenon, known as confinement, cannot be derived by the 
perturbative analysis. 

From experimental data of particles one can see three lightest hadronic 
states 7T°, 7r + , 7T~ belonging to the triplet of the SU(2) multiplet of u, d 
flavors and 7r°, 7r+, 7r~, rj, K + , K~ , K°, K° belonging to the octet of the 
577(3) multiplet of u, d, s flavors with the masses below the mass of the 
higher lying p meson. At present it is generally assumed to be a result of the 
spontaneous chiral symmetry breaking which occurs in the massless QCD 
Lagrangian. The non-vanishing quark masses break the chiral symmetry of 
the Lagrangian explicitly which gives rise to non-zero masses of the quasi- 
Goldstone bosons. The particles in the multiplets are identified with these 
quasi-Goldstone bosons. Their masses are proportional to the square root of 
the quark masses. This phenomenon, known as the chiral symmetry breaking, 
is a truly non-perturbative effect. The chiral symmetry that is broken at low 
energies is, however, expected to be restored if the energies are increased. This 
scenario corresponds to a phase transition. An order parameter indicating 
this phase transition is the chiral condensate (ipip), which is zero in the 
chirally symmetric phase and picks up a non-zero value once we pass to the 
chirally broken phase. 

Lattice QCD is a theory that provides a non-perturbative approach from 
the first principles. One introduces a space-time lattice and formulates the 
theory by a lattice regularized path integral. Quarks and gluons become regu- 
larized fields. Quarks are represented by Grassmann spinor fields attached to 
the lattice sites and the gauge fields are represented by parallel transporters 
defined on the lattice links. The regularization introduces a momentum cut- 
off n /a, where a is the lattice spacing. This represents a statistical model 
that is subject to the Monte Carlo simulations on a computer. The contin- 
uum limit is recovered if one sends the lattice spacing to and the number 
of lattice sites to infinity, keeping the physical volume fixed. 

When addressing the phenomenon of the chiral symmetry breaking and 
the physics of the light quarks, it is vital to have a discretization of fermions 
which preserves the chiral symmetry at finite lattice spacing. In this thesis 
we will describe results obtained with various formulations of such chiral 
fermions. At present all these formulations are tedious to simulate in QCD 
due to their computational overhead compared to the non-chiral discretiza- 
tions such as the Wilson fermions. Therefore all simulations are carried out at 
quite small physical volumes and the values of the meson masses higher than 
one would expect in nature. The latter is reflected for instance in the pion 
mass. Hence the simulated data have to be necessarily extrapolated to the 
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small physical quark masses. In addition the quark loops are also neglected 
which leads to quite difficult accountable systematic uncertainties. A way to 
extract physical information from the QCD simulations is to use theoretical 
predictions provided by a low energy effective theory. In the context of QCD 
this is the Chiral Perturbation Theory (xPT). 

Owing to the mass gap between the quasi- Goldstone bosons and the first 
heavier mesons (the p-meson) one can integrate out all higher excitations of 
the mass spectrum of QCD leaving the theory defined solely in terms of the 
quasi-Goldstone bosons. The resulting theory is a low energy effective theory, 
which is called the %PT. 

An attractive feature for the numerical simulations is that %PT can be for- 
mulated in a finite volume and some analytical predictions can be computed 
which are the same as in the infinite volume. Moreover since the present lat- 
tice simulations are done without quark loops this can also be accounted for 
in the quenched which provides an analytical guide-line for quenched 
lattice simulations. In particular the analytical predictions in the framework 
of xPT have been worked out for the regime of the soft pion where the Comp- 
ton wavelength cannot be accommodated in the finite volume. This regime, 
though representing an unphysical situation, provides a powerful tool for the 
lattice simulations with chiral fermions. The main advantage is that the finite 
size effects can be calculated in this setting using the ;\PT an d the low energy 
constants appearing in the corresponding analytical expressions are those of 
the infinite volume. The control over the finite size effects, which one other- 
wise wishes to suppress, provides us with physical information. Since in this 
regime one does not need to use large lattices, it is viewed to be promising 
for obtaining the low energy constants parameterizing the ;\PT Lagrangian. 
This can be achieved by using the formalism of the chiral fermions. The use 
of the chiral fermions allows us to perform calculations at small meson masses 
and thereby approach their physical values, i.e. the so-called chiral limit of 
observables of interest. 

In Chapter 2 we give an introduction and fix our notations in QCD, xPT 
and chiral Random Matrix Theory. Lattice QCD will also be briefly dis- 
cussed. For completeness the traditional formulations of lattice fermions will 
be reviewed. We discuss various formulations of the chiral lattice fermions 
realized through the Ginsparg- Wilson relation. In particular we consider the 
conventional formulation of the overlap fermions and an improved formula- 
tion based on the perfect action — the hypercube fermions. The Chapter 3 is 
concerned with simulation aspects of the Ginsparg- Wilson fermions and also 
quenched gauge simulations are briefly discussed. The Chapter 4 deals with 
the numerical aspects of simulation of the hypercube fermions. In Chapter 5 
we present our results for the dependence of the pion mass versus the bare 
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quark mass. In Chapter 6 the results for the distributions of the eigenvalues 
of the overlap Dirac operator are considered in conjunction with the chiral 
Random Matrix Theory. In Chapter 7 the results for the topological suscep- 
tibility, computed with the chiral fermions, are addressed. In Chapter 8 the 
two-point mesonic functions are addressed. We discuss the feasibility to ex- 
tract the physical information from the axial and pseudo-scalar correlators. 
In Chapter 9 we present a study of a gauge action which may be useful for 
the future lattice simulations in the e-regime. In Chapter 10 we draw the 
conclusions and give an outlook. 
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Chapter 2 

Theoretical background 



2.1 QCD at first glance. Gluons and quarks 

In this Section we introduce the basic notation of QCD. Its degrees of freedom 
are carried by quarks and gluons. Gluons are described by the non-Abelian 
SU(3) gauge group while quarks belong to its fundamental representation. 
We will consider QCD in Euclidean space that will be the relevant theory for 
the results discussed in the thesis. The Euclidean formulation is especially 
useful for the lattice simulations since the theory becomes a statistical en- 
semble with the Boltzmann weight exp (— Sqcd), where Sqcd is the QCD 
action. The use of the Euclidean formulation is justified since there exists an 
analytical continuation of the Euclidean Green functions to the Minkowski 
space [3, lif • 

The quarks are represented by spinor fields in the fundamental representa- 
tion of the SU (3) color group, and the gluon vector fields (// = 1, . . . , 4) 
take their values in the Lie algebra of the SU (3) color group. The quarks 
come in six flavors ipf, where / = 1, . . . , 6. We decompose the matrix valued 
vector potential in fields A®, where a = 1, . . . , N% — 1 with N c being the 
number of colors, 



g is the bare coupling constant and the matrices T a are the generators of the 
gauge group SU(3), which are traceless and Hermitian. They are normalized 




(2.1) 



a=l 



by 



Ti(T a T b ) = -8 ab . 
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The generators satisfy the commutation relation 

T'T 



z/abcy , (2.2) 



where f a b c are structure constants of the SU(3) Lie group which can be 
chosen to be real. The values of the structure constants for the special choice 
of generators as the Gell-Mann matrices are given in Ref. 

We split the QCD action into two parts, Sqcd = Sym + SV> with the 
action of quarks given by 

S F = V / d 4 xip f (V + m f )iJ)f . (2.3) 
f=1 Jv 

Nf is the number of quark flavors, SVm is the Euclidean Yang-Mills action, 
T> is the massless Dirac operator and m/ is the quark mass corresponding 
to the flavor /. The integration is done over the Euclidean four-volume V = 
L 3 x T. The boundary conditions are taken to be periodic up to a gauge 
transformation for the gluons, and antiperiodic for the quarks. This defines 
the theory on a four dimensional torus. The QCD Dirac operator is given by 

V = 7^ + A„) . (2.4) 

The 7 M are the Euclidean Dirac matrices with {7^, 7^} = 25^. We chose them 
such that the Dirac operator is anti-Hermitian, X>t = —T>. In particular we 
use the chiral representation in which 75 = 71727374 = diag (1, 1, —1, — 1) 
(see for instance Ref. 0). 

The Euclidean Yang-Mills gauge action is given by 

Sym = \ [ d*xF a F« . (2.5) 



4 



V 



In this formula, and throughout this thesis, the summation over the repeated 
indices is assumed if not stated otherwise. The field strength tensor F^ v 
carries one color a index running from 1 . . . AT? — 1 and two Lorentz indices 
H,u = 1...4, 



F; v = d,A a u - d v Al + gU^Al . (2.6) 

The theory is finally defined by its partition function in the form of the 
functional integral over the fermion fields ipf and ipf and the gauge vector 
potential A^ in the Euclidean space, 

Nf 

Z QCD (m f , 9) = I Y[[DA^ ]\[Di; f Di; f } e s Q cu^ f ^ f ,A^[A] _ ^ 
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The parameter 6 is a real vacuum angle and the topological charge v is 
defined by 

v\A\ = -^J^d 4 x e a/3pa Ti(F a(3 F p(T ) (2.8) 

where e a p prj is the antisymmetric unit tensor of rank four. v[A] is the winding 
number of a field configuration and takes integer values for finite action field 
configurations. The field ipf is a four component spinor field and ipf is a 
conjugate spinor field. Integrating out the fermion degrees of freedom we 
arrive at 

Z QCD (m f ,9) = I Y[DA„ Y[det{V + m f )e- s ™W +i9 »W. (2.9) 

The integral over all field configurations J f| [-DA M ] includes a sum over all 
topological sectors v. 

The partition function in the sector of topological charge v is obtained 
by the Fourier transform 

Z? CD (m f ) = — r dQe" ive Z^[m j} Q) . (2.10) 

2.2 Chiral symmetry breaking and eigenval- 
ues of the Dirac operator 

2.2.1 Implications of the chiral symmetry for the 
structure of the Dirac operator 

In this Section we discuss the chiral symmetry in QCD and its implication 
for the spectrum of the Dirac operator. 

We start from the definitions of right-handed and left-handed spinors 
ip R / L which are defined as projections of ip and ijj 

P (+M = ^(1±75), 

1 = p(+)+p(-), p(+)p(-) = 0, (2.11) 

where the spinor describes Nf Dirac spinor species. 

The massless Dirac operator satisfies the anticommutation relation 

{ l5 ,V} = 0. (2.12) 
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Using the decomposition of the fermion fields if), if) 

if) = if) R + if> L , ^ = ^ + ^ L (2.13) 
we obtain for the fermion action 

S F = J d 4 x(ip R Vif) R + i) L Vip L + ip R Mip L + i) L Mif) R ^j , (2.14) 

where M = diag(m 1; . . . , rrif). 

Relation ()2.12|) can be seen also as an expression of the chiral symmetry 
of the fermion action for every fermion flavor if>f, 

if) f -> e iai5 ip f , if) f -»■ %f) f e ia ^ , (2.15) 

where a is a real angle of the chiral rotation. The quark action ()2.14j) remains 
invariant under a chiral rotation if M = 0. 

One can write down an eigenvalue equation for T> as 

V[A]if) n = i\ n if) n , X n eR, (2.16) 

where the eigenvalues and eigenfunctions depend on the gauge field. Using 
Eq. (|2.12j) one can show that the nonzero eigenvalues of T> occur in pairs 
±i\ n with eigenfunctions if> n and ■j5ip n - Indeed, if i\ n is an eigenvalue of V 
and if) n is the corresponding eigenfunction, then Eq. (|2.12|) implies 

£>75^n = -Ifplfln = -iKlS^n , (2.17) 

hence is a ls° an eigenfunction of T> with the eigenvalue —i\ n . There 

can also be eigenvalues equal to zero, A n = 0. Because of Eq. ()2.12|) the 
corresponding eigenfunctions can be arranged to be simultaneous eigenfunc- 
tions of V and 75 with eigenvalues ±1 of 75, i.e. these states have a definite 
chirality. Denoting the number of zero eigenvalues per flavor / with positive 
and negative chirality by iV + and AL, respectively, the Atiyah-Singer index 
theorem states that 

v = N f (N+ - N-) (2.18) 

is a topological invariant that does not change under continuous deformations 
of the gauge field. Note that the Atiyah-Singer index theorem identifies the 
definition of the index given by gauge fields (|2.8j) and the definition given by 
the index of the Dirac operator, i.e. the fermionic index. 

In a chiral basis with jsipn = ±if>n^ L , one can use Eq. ()2.12j) to show 
that (ip R \V\il)^) = = (V^l^l^n) for all m and n. From this property 
and the fact that D is anti-Hermitian, it follows that the projection of the 
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Dirac operator onto the subspace spanned by the zero modes has the matrix 
structure 

* m - ( T ) • ^ 

This off-diagonal block structure is characteristic for systems with chiral 
symmetry. The matrix W has dimension NfN + x NfN_, and the matrix T> 
in Eq. (|2.19|) has \u\ eigenvalues equal to zero. 

Using formula ()2.19jl we can write down the contribution of the zero 
modes to the fermion action as 

Sf = J d 4 x(^p R iW j tp R + ^p L iWip L + i) R Mip L + i) L Mi) R ^j . (2.20) 

In the chiral limit where all rrtf = 0, the fermion actions ()2.14j) as well 
as ()2.2()j) are invariant under the global transformations 

i) L -> Lib L , %b L -> i) L L~ l , 

(2.21) 

L G U L (N f ) and i? G U R (N f ), where U L (N f ) and U R (N f ) are independent 
unitary groups. If the number of right-handed states = NfN + is equal 
to the number of left-handed states Nl = NfN_, the symmetry is thus 
XJji(Nf) ® XJ L (Nf). However, if N R ^ N L , i.e. if u ^ 0, the axial symmetry 
group Ua(1), which is a U(l) subgroup of XJ R (Nj) <8> XJ L (Nf), is broken by 
instantons. Invariance under a second subgroup U(l) of XJ R (Nf) <S> Ul(./V/), 
the group Uy(l), corresponds to the conservation of the baryon number. 
The full flavor symmetry group in the chiral limit is thus given by G = 
SV R (N f )®SV L (N f ). 

Since G is a symmetry group for the massless fermion action one can 
write down the Noether currents 



= VTu^r^T^. (2.22) 



(1 + 75) 



T a is a generator of the SU {Nf) group, where a — 1, ... , Nj — 1. One obtains 
conserved Noether charges by the usual procedure, 

Q R J L = Jd 3 xjS L (x). (2.23) 
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We construct generators of the axial SUa(-A^) and vector subgroup SXJy(Nf) 
by linear combination of the Noether currents J L and J R , 

jA jR _ jL 

J£ = Ja, + Ja,- (2-24) 

The other two currents which are related to the Uy(l) and Ua(1) symmetries 
are given by 

= ^757/,^ , 
Jl = (2-25) 

As we mentioned before, for v 7^ the axial current is not conserved. Its 
divergence is given by 

Tr[F^(x)F pCT (x)] . (2.26) 

One can also see that the axial SU (Nf) subgroup with L = R~ l is broken 
explicitly by the mass term. On the other hand the SU(A^) vector subgroup 
(with L = R) is not broken for degenerate quark masses (rrif = m for all /), 
but breaks explicitly for different quark masses (nif 7^ rxif for some / 7^ /'). 

If the quark masses were zero, the action ()2.20|) would have an exact axial 
flavor SU(iV/) <g> SXJ(Nf) symmetry. However, there is strong evidence that 
the axial SU(Nf) subgroup with L = R~ x is broken spontaneously with the 
following spontaneous symmetry breaking pattern of the chiral group, 

SV R (N f ) <g> SV L (N f ) -> SV v (Nf) . (2.27) 

This implies that the ground state of the quantum system being invariant 
with respect to the SUy(iV/) is, however, not invariant with respect to the 
SU^iVy) symmetry. According to the Goldstone theorem spontaneous break- 
ing of a Lie group G down to a Lie group H gives rise to Goldstone boson 
fields 0,0. The Goldstone manifold is a coset space G/H = SXJ(Nf), and so 
there are Nj — 1 Goldstone bosons. 

The fact that the pattern of the spontaneous symmetry breaking has the 
form S\J R (Nf) ® SU^iV/) to SUy(iV/) is supported firstly by the absence 
of the parity doublets in the hadron spectrum. Indeed in the mass spec- 
trum one can observe three light pseudo-scalar particles, the pions (which 
are the triplet of an approximate SU(2) meson multiplet) as well as some- 
what heavier pseudo-scalars — the four kaons K + , K°, K°, K~ and the n 
-meson (belonging to the approximate SU(3) meson multiplet). Along with 
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the eight approximate Goldstone bosons one observes a heavier rj meson. Its 
mass is too large to be considered as a result of the spontaneous symmetry- 
breaking of the axial Ua(1) subgroup. According to 't Hooft's idea its mass 
receives contributions from instantons due to a relation of the derivative of 
the anomaly axial current to the topological charge [1, 

Secondly the value of the chiral condensate (ipip) = (ip R ip L ) + (ip L ip R ) 
would be zero if the ground state was axial-flavor symmetric. However, lattice 
QCD simulations indicate that (ipip) is non-zero. The quantity (tpip) is only 
invariant if L = R, i.e. the vacuum state is symmetric under the flavor group 
H = SUy(Nf). Thus, the vector symmetries are unbroken, whereas the axial 
symmetries are maximally broken. If we now add to this scenario of the 
spontaneous chiral symmetry breaking the explicit symmetry breaking by 
the mass terms of light quarks, then we arrive at a consistent picture for the 
observable approximate SU{2) triplet and 577(3) octet. Comparing masses 
of the pions in the triplet and pions, kaons and n meson in the octet, we see 
that the explicit symmetry breaking is larger for the octet. This implies that 
the mass of the strange quark in the QCD Lagrangian is large compared to 
masses of the u and d quarks, which one often assumes to be degenerate. 



2.2.2 Spectral density of the Dirac operator 

In terms of the eigenvalues of the Dirac operator, the QCD partition func- 
tion f)2.9j) can be rewritten as 

+00 N f 

Z QCB (m f ,6) = e ^n m /' / DA^\{{\l + m))e- s ™W ; ( 2 .28) 

u=-oo f=l ^ u k>l 

where f DA^ denotes the path integral over field configurations with topo- 
logical charge v. The eigenvalues of the Dirac operator, X k , are introduced in 
Eq. (jHSD- 

The spectral density of the Dirac operator is given by 

KA) = (J>(A-A n )) , (2.29) 

n 

where the expectation value is over gauge fields weighted by the full QCD 
action. The spectral density is important because of its relation to the order 
parameter for spontaneous chiral symmetry breaking, the chiral condensate 
£ = -$V>- ' _ 
It was shown by Banks and Casher jjJJ] that 



S = lim lim lim np(e)/V . (2.30) 

e— >0 m^O V-^oo 
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It is important that the limits are taken in the order indicated. (In the nor- 
malization of Eq. (|2.29j) the spectral density is proportional to the volume, 
so the explicit factor of 1/V in ()2.3()j) is compensated to yield a finite result.) 
This relation can be readily derived. The chiral condensate is given by 

(W) = - lim lim -l-d m \ogZ QCD (m,6) . (2.31) 

m— »0 V-*oo V is f 

If we now apply this definition to the partition function ()2.28|) then we arrive 
at 

(W) = - lim lim ( - V — > , (2.32) 

n>l n 

where the contribution of the zero modes, \i/\/(rnV), was dropped assuming 
that it becomes asymptotically small in the V — > oo limit. Now taking the 
chiral limit we obtain the Banks-Casher relation ()2.30j) . 



2.3 xPT as a low energy effective theory of 
QCD 

Due to spontaneous chiral symmetry breaking, at low energies the QCD 
dynamics is dominated by Goldstone bosons p, u| the pions for Nf = 2. 
It is possible to use a low-ener gy effective description that only involves 
the quasi- Goldstone boson fields ll lij]. Chiral perturbation theory 14, 15| 
provides a systematic low-energy expansion that predicts the pion dynamics 
based on symmetry principles and a number of low-energy parameters (like 
the pion decay constant F n and the chiral condensate S) whose values can 
be determined either from experiments or from lattice QCD calculations. 



2.3.1 The chiral Lagrangian and its low energy con- 
stants 

The Lagrangian of the low-energy effective theory is determined by the global 
symmetries of the underlying theory. Considering the low energy effective 
theories one has to think of an energy scale hierarchy and the corresponding 
theory which is relevant for a given energy interval. As we reduce the energy 
scale, some degrees of freedom, which were important for higher energies, 
may become irrelevant for the lower energies. 

The low energy effective Lagrangian is constructed as a linear combina- 
tion of all terms which are local, Lorentz invariant and compatible with the 
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flavor symmetry. The constants in such an expansion parameterize our lack 
of knowledge about the underlying theory. 

If we write down a low energy effective theory for QCD then the corre- 
sponding flavor symmetry includes the chiral symmetry. It is broken down 
to SU (Nj) and therefore the Goldstone bosons belong to the coset space 
SUfl(JV)) <g> SV L (N f )/SV v (Nf), which is again SU(iV». If we are interested 
only in the lightest excitations then the Lagrangian will contain only Gold- 
stone fields U(x). Under global chiral rotations they transform as 

U(x)' = LU{x)I$ , (2.33) 

where L and R axe matrices of the chiral transformation belonging to 
SU L (Nf) and SU R (Nf), respectively. To fully determine the theory we need 
to write down its partition function and specify the measure. The partition 
function of the chiral perturbation theory takes the following general form 

Z xPT = J [DU{x)\ exp (- J d 4 xC[U]^j , (2.34) 

where [DU(x)\ represents the Haar measure and the integration is done over 
the flavor group SU(Nf) in each space point x. The Haar measure is a left- 
arid right-invariant measure, i.e. 

/ dU f(ttU) = [ dU f(UQ) = [ dU f(U) , (2.35) 

JSU(Nf) JSU(Nf) JSU(Nf) 

for any function f(U) and for any matrix f2 G SU(Nf). It is convenient to 
normalize the measure such that 

/ dU =1 . (2.36) 

JSU(Nf) 

The interaction among the Goldstone bosons is weak at low energies and 
decreases as the energy is reduced. For the effective Lagrangian it implies that 
it can be represented as an expansion in derivatives. On the other hand the 
quark masses can be viewed as parameters which characterize how strongly 
the chiral symmetry is broken. 

To a good absolute accuracy the up and down quarks can be regarded to 
have nearly the same mass. So for Nf = 2 we can choose the mass matrix 
M. = ml. One can write a systematic low-energy expansion in powers of the 
pion momentum p and mass tti/^o — 2,mYj / over the cut-off of the effective 
theory, Aqcd — 4:7rF n . In the standard chiral expansion both quantities, A ^ g 
and j^° D , are taken to be of the same order. Because of these counting rules 
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the chiral expansion cannot be considered as a Taylor expansion. Here we 
give the result for the leading term of the pion effective Lagrangian in the 9 
vacuum, 

£(2) = S-Trid^d^U) - ^Tr(e™tMU* + Ue~^ ' M ] ) . (2.37) 

The first term on the right-hand side is chirally invariant. Its prefactor is given 
by the constant F n related to the pion decay constant which determines the 
strength of the interaction between the Goldstone bosons. The second term 
is the chiral symmetry breaking mass term which contains the quark mass 
matrix. Under chiral transformations ()2.33|) this term transforms as 

Tr {e^ MU' 1 + U'e' l ^f M r ) = Tr{e^iMBrftf + LUR ] e~ i '^ M ] ) . (2.38) 

Since we assume the two quarks to have identical masses, the Lagrangian 
is invariant under SU(2) flavor rotations for which R = L. For a general 
diagonal mass matrix with Nf flavors the flavor symmetry is reduced to 

11/ ',/•(!)/• 

The constants F n and £ determine the low-energy dynamics at leading 
order and enter the effective theory as free parameters. 
It can be shown that in the chiral limit 

(0\J*\Tr + (p))=tV2F 7T p fl . (2.39) 

The same formula at physical values of quark masses is in fact the definition 
of the pion decay constant F^ ecay . So the physical value of F^ ccay being F n 
in the leading order is modified by 0(mj) ~ 0(p A ) corrections in the next- 
to-leading order, 

Ff cay = F 7r + 0{m 2 f ,p 4 ) . (2.40) 

The pion decay constant is experimentally known from the process 
^ —>■ ^ + u and it is approximately 93 MeV. If we formulate the chiral 
Lagrangian for the 577(3) flavor group — i.e. including the s-quark — then 
in addition to the pion decay constant we will have also the kaon decay con- 
stant. In principal one could relate it to the prefactor of the kinetic term 
of the chiral Lagrangian. The kaon decay constant is about 114 MeV. The 
difference between the pion decay constant and the kaon decay constant is a 
typical 0(p 4 ) effect, which is beyond the lowest order. However, the effects of 
order C(p 4 ) are expected to be larger in the case of the kaon decay constant. 
Therefore the most natural determination of the prefactor at 0(p 2 ) is still 
provided by the pion decay constant. 
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Contrary to the decay constants, the quark condensate is not directly 
related to any physical observable. It is the product m q B, where B = T,/F%, 
that can be constrained by means of experimental data. To the leading order 
we have 

mlo = (m u + m d )B, 

m 2 K o = {m d + m s )B . (2.41) 

The analogous equation for m 2 s , where n$ denotes the octet component of the 
rj meson, contains no new parameters and gives rise to a consistency relation, 

3m 2 vs = Am 2 K0 - m' . (2.42) 

This is the famous Gell-Mann-Okubo mass formula, well satisfied by ex- 
perimental data under the assumption m v = m Va . The validity of the Gell- 
Mann-Okubo formula provides a significant a posteriori check that the 0(m 2 ) 
corrections to Eq. ()2.41j) are small. 

Up to these two low-energy constants the Goldstone boson dynamics at 
the leading order of the chiral perturbation theory is completely determined 
by chiral symmetry. At higher energies additional terms arise in the effective 
theory. Again, they are restricted by chiral symmetry and they contain new 
low energy constants — the Gasser-Leutwyler coefficients. 



2.3.2 The p— and e— expansion of xPT 

Let the system now to be placed in a four dimensional box of length L. We 
assume from now on that the volume is large with respect to the QCD scale 
set by the 1/Aqcd- More precisely, the allowed momenta should not take 
values higher than the scale of the chiral symmetry breaking in the infinite 
physical volume given by Aqcd = ^F n , 

p ~ — < 4vrF^ . (2.43) 

As a general property, no chiral symmetry breaking occurs in a finite volume. 
However, this condition separates the hard momenta from the soft ones and 
thereby provides the corresponding framework for calculations in the finite 
volume. 

The inverse product l/(F n L) is not the only relevant parameter. Also the 
relative size of and 1/L is important, as it was discussed in the previous 
Section. Indeed the %PT considers an expansion in two parameters: the mass 
of the quasi-Goldstone boson m n and its momentum p over the scale 47iF n . 
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In the finite volume the corresponding momenta are quantized in units of 
2ii/L. The Euclidean partition function for the xPT, 

Z{V,M) = J[dU}exp(^-J d^xCp^j , (2.44) 

involves three energy scales: Aqcd, -M and 1/L or p = l^njL. 

In the chiral expansion, U = exp(iy/2r)(x) / F n ) is expanded around the 
classical solution 

U(x) = 1 + iV2n(x)/F^ + ... . (2.45) 

The field rj(x) belongs to the Lie algebra of SU(3) and it describes fields of 
the pseudo-scalar mesons. 

Inserting Eq. ()2.45|) in the expression for the Lagrangian ()2.37|) and inte- 
grating over the volume we arrive at the action 

/ d 4 x £ (2) = —NfF^mlV + lf d 4 x Tr {{d^f + m^ 2 ) . (2.46) 
Jv ^ ^ Jv 

Since rj(x) is defined on the torus V = L 4 we may expand it in terms of 
periodic plane waves 

oo oo 

v( x ) = z2 " zJ q n u n (x) + h.c. , n = (ni, . . . ,n 4 ) , 
m=— oo ?i4=— oo 

u n (x) = exp(2ni nx / L) , (2.47) 

where h.c. represents the Hermitian conjugate of the preceding term. Insert- 
ing this expression in Eq. 1)2.46)1 we see that the zero mode degrees of freedom 
go of course do not appear in the kinetic term. This fact manifests itself also 
in the expression for the pion propagator which develops a pole, 

G(l _ j/) = "°W".W + ... (2 .48) 
7774 V 

So now one can distinguish several cases. If the linear size L is much 
larger than the Compton wavelength of the pions, l/m n , then the system 
hardly feels the finite volume and we have the usual chiral expansion. As the 
size of the box or the pion mass becomes smaller, finite size effects become 
significant, but provided that m n L > 1, ordinary perturbation theory is still 
applicable. This is the so-called p-expansion [la] with the counting rules 

\q n \ m w p 1 

r^j r*j r*j 

F-k A.QCD A.QCD LFtt 
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If the chiral limit is approached further in such a way that the Compton 
wavelength of the pion is larger than the box size (L < l/m n ), the conven- 
tional p-expansion eventually breaks down due to the propagation of pions 
with zero momenta . Indeed, according to Eq. ()2.48|) the pion propagator 



contains a pole. Therefore in order to arrive at a representation which would 
be valid in this regime of small quark masses (or lengths L) one needs to 
reorder the chiral perturbation series by summing graphs which involve an 
arbitrary number of zero mode propagators. 

An approp riate expansion for this regime is the so-called e- expansion, in 



An appro] 
which FuJil 



Aqcd A 2 QC d 

Note that S and F n , as well as other low energy constants, appearing in this 
counting scheme are those of the infinite volume. Using Eq. ()2.43|) the above 
relation implies for the e-regime 

vrw < j < 4-71-^ . 

The difficulty in this regime comes from the fact that the zero modes have 
to be accounted for exactly, while the chiral perturbation theory applies to 
the non-zero mode integration, since ~ 0(e 2 ) for \p\ ^ 0. One treats the 
zero modes of the pion as a collective coordinate by factorizing U(x) into a 
constant collective field U and the pion fluctuations £(x), 

U(x) = tfoexpz^Mj . (2.49) 

The field £ (x) can be represented in terms of the plane waves with momenta 
\p\ larger than zero. Uo is a SU(Nf) matrix independent of x. 
To leading order, the partition function reads 



2(9, m q ) 



d£dU exv\ I ^xTr (d^(x)d^(x)) x (2.50) 

SU{N f ) 2 Jy 

exp — Tr (Me ie/N fU + h.c.) , 
M = diag(mi, . . . , m Nf ) . 

It is also interesting to consider expectation values in sectors of fixed 
topology [2^|. Fourier-transforming in 9 and dropping the non-zero mode 



17 



Chapter 2 Theoretical background 



contributions, we obtain 

Z v ( mi ,...,m Nf ) = [ dUoidetUoY exp(^Tr( MU + U f Q MA) 

JsU(N f ) \ 1 V V 

or for the degenerate quark masses with M. = diag(m g , . . . , m q ) 

Z v (m q ) = [ dU (detU o y exp ( H^Tr (u + U}S) . (2.51) 
Jsu(N f ) V 1 v > 7 

If F^L ^> 1, the p- and e-expansions should match in the range of quark 
masses such that m n L ~ 1. In this regime m q YM ~ (F n L) 2 (m n L) 2 ^> 1, so 
the results in the e-expansion should reproduce those of the p-expansion in 
the limit of large m q Y,V . 

2.3.3 Quenched xPT: first order expressions for the 
axial 

-vector correlation function 

We start by writing down the expression for the axial- vector current, 

A^t, x) = $(t, x)^T^(t, x) . (2.52) 

The two-point correlation function is as usually given by the expectation 
value 



{A^x)A u {y)) - Z ^(m q ,9) ^ j 

or in terms of the propagator (V + m g )~ 1 (x, y), 

(Al(x)A b u (y)) = 

-Tr(T a T b )(Tr[ l5 ^(V + m q )-\x,y)^(V + m q )- 1 (y,x)\) . 1 ' ) 

This expression can be obtained using the partition function for QCD in the 
presence of auxiliary external fields a M (x), the so-called generating functional. 
For the sake of the following discussion we will also add auxiliary fields for 
the pseudo-scalar density P(x) = ipi^ip , 



Z(a fl ,P)=J DA^D^Dipexp -Sqcd + J ^^(7^75^ + hbP)^ 



2.55) 
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If we compute the second functional derivative of the partition function 
Zfan, P) with respect to the axial current a M and put a M and P to zero 
afterwards, we recover the expression (|2.5H|) for the axial current correlation 
function. This new partition function ()2.55|) is now subject to its low energy 
approximation, i.e. one integrates out the quark degrees of freedom leaving 
only the mesons. The leading order expression of the chiral Lagrangian of 
the low energy effective theory in the presence of a M and P takes the form 



C^[U,a„P] 



(2.56) 



2 



f(M + iP{x))U ] + Ue N f(M ] - iP\x)) 



where 



V jJJ = dfJJ — ia^U — iUa^ , M. — diag(m g , . . . , m q ) . 

We see that the partial derivatives in the expression (|2.H7jl were promoted 
to some kind of covariant derivatives. Note that P(x) is assumed to be of the 
same order as m q in the chiral expansion, whereas a^(x) is counted as the 
pion momenta p . 

Now from this expression one can construct the partition function. Taking 
then two functional derivatives with respect to the axial current a M (x) or 
two functional derivatives with respect to the pseudo-scalar density P(x) we 
obtain the expression for the axial correlation function and pseudo-scalar 
correlation function in %PT. As far as full QCD is concerned, all two-point 
functions had been calculated before for the e-regime in Ref. [20]. 

The QCD partition function (|2.9jl includes a fermion determinant. It 
accounts for the creation and annihilation of the quark- ant iquark pairs. It 
turned out that in computer simulations the calculation of the determinant 
slows down the updating of the gauge fields tremendously since one has to 
recompute it for every gauge configuration. Therefore a crude approximation 
to QCD that goes under the name of quenched approximation is used in part 
of the contemporary simulations. It neglects the fermion determinant in the 
QCD partition function ()2.9|) and in the expectation values of the observables 
calculated with it. Despite the crudeness of the approximation it is useful to 
employ it as an estimate for the observables of interest, keeping in mind 
the systematic uncertainties of about 10 — 15% 24, 2^. The latter depends 



of course on the observable of interest. We aim, however, at the quenched 
formulation of the chiral perturbation theory, i.e. a low energy effective the- 
ory where the quark loops were eliminated on the level of QCD. In order 
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to get rid of the quark loops one can add quarks with bosonic statistics to 
the fermion part of the QCD Lagrangian (|2.H|) keeping the external currents 
coupled to the quarks with the proper statistics. The latter are called the 
"valence" quarks. This results in the cancellation of the diagrams with inter- 
nal loops. The corresponding partition function with the sources Jj coupled 
to quark bilinears is, 

...,JnJ-J [DA,] det{v + mq)Nv e 

(2.57) 

The 2N V quarks (N v with usual fermionic statistics and N v additional with 
bosonic statistics) lie in the fundamental representation of the graded group 
SU(N V \N V ). For massless quarks, the corresponding full QCD Lagrangian 
exhibits a graded symmetry SU(N V \N V ) L ® SU(N V \N V ) R <g> U(l) v that is 
assumed to be spontaneously broken down to SU(N V \N V ) V ® U(l) y. This is 
the so-called "supersymmetric" formulation of quenched QCD [26(. Using it 
as a starting point one obtains a low energy theory without quark loops. 

The "supersymmetric" chiral Lagrangian in the leading order is given by 

4usy = X Str i d v U ~% U ) - ^Str(W + U-'Ug 1 ) + (2.58) 



9 



where Str denotes the "supertrace" , $o = ^Str[— i log(Z7)] and 

U 9 = exp(i9/N v )I Nv +I Nv . (2.59) 

Here, In v is the identity matrix in the fermion-fermion block of "physical" 
Goldstone bosons and zero otherwise, while Im v is the identity in the boson- 
boson block and zero elsewhere. Apart from the familiar low energy constants 
F n and S — which now, however, take their quenched values — we have 
also two new parameters m and a, which are artifacts of the quenched 
approximation. 

We note that there is also another method to formulate the quenched 
xPT — the replica method j2?| . Here one considers valence quarks to be 
part of iV quarks in total. Considering now N as a parameter and sending 
it to we obtain the replica limit. The equivalence of the two methods was 
verified at the perturbative level of the chiral Lagrangian |27i |. 

An important difference of quenched from the full theory is the 

unavoidable presence of a singlet scalar field $ 0) which cannot be decoupled 
in this case. The singularities of the propagator of the $o held are responsible 
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for the quenching artifacts. Also the formulation of the counting rules for the 
quenched becomes more involved. These issues, as well as the two-point 
functions for the e- regime in the quenched approximation, were addressed 
in Refs. 0Hj. 



For the purpose of this thesis we consider the expression for the axial 
correlation function with one flavor. We define the bare axial current with 
one flavor at momentum p = and Euclidean time t, 

Mt) = ^TsT^fr x) • (2.60) 

X 

The formula for the axial correlation function in a volume L 3 x T, to the first 
order in quenched xPT, reads ji^ 1 



2-f^ + 2m q E W] (z q )T-h 1 (r) S j 



hi(r) = -\t 2 -t + -), r = 4 , 
2 V 6/ ' T 

K{z g ) = X (z q [I v (z q )K v (z q ) + I v+x (z q )K v _ x (z q )] + j")(2.61) 

7„ and K v are modified Bessel functions, z q = rrigT.V, v is again the topo- 
logical charge, and h\ is a purely kinematic function. The latter shows that 
this correlation function is given by a parabola in t with its minimum at 
T/2. This behavior is qualitatively different from the infinite volume, where 
the correlations decay exponentially (resp. as a cosh function in a periodic 
volume). However, the small and the infinite volume have in common that 
in both cases the quenched axial correlation function to the first order only 
depends on the low energy constants and S. 

2.3.4 Zero mode contributions to the pseudo-scalar 
correlation function 

We start off by writing the expressions of the pseudo-scalar densities 
V°(x),V a (x) 

p° = ipi^ip , 

V a = i>T a i lb ^, a = l,...,N 2 f -l, (2.62) 



1 We have an extra factor of 2 compared to Ref. [2j| because this definition of the 
current is considered only for one flavor. 
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where T a are generators of the group SU(Nf) and a runs from 1 to Nj — 1. 
Let I, J run from to Nj — 1 and T° = 1, then the correlation function 
(V 1 (x)V J (y))u takes the general form 

(P / (x)P J d/)), = -Tr(T / T J )P 1 (x,i/) + Tr(r / )Tr(r J )P 2 (x, 2 /) . (2.63) 

V\ and Vi are the connected and disconnected parts of the correlation func- 
tion, 

Vi(x,y) = Tr[i^ 5 (V + m q )~ 1 (x,y)i-f 5 (V + m q )~ 1 (y, x)] , 

V 2 {x,y) = Tr[i l5 {V + m q )- l {x,x))Tr[i l5 {V + m q )- 1 {y,y)) , (2.64) 

where P _1 (x,y) is the propagator of the Dirac operator. The correlation 
function can be computed from the generating functional which we discussed 
in the previous Subsection. Now one has to calculate two derivatives of the 
generating functional with respect to the auxiliary pseudo-scalar field P(x) 
introduced in Eq. ()2.55|) and put all the auxiliary fields to zero at the end. 
Applying the same technique to the generating functional of the low energy 
effective theory we obtain the corresponding pseudo-scalar correlation func- 
tion in the quenched xPT (q%PT). 

We consider the pseudo-scalar correlation function for the definite topo- 
logical sector v. Application of the spectral representation of the quark 
propagator to the pseudo-scalar correlation function was first considered in 
Ref . [3i| . By employing the spectral representation of the quark propagator, 
it is clear that the correlator (V 1 (x)V J \y)) v contains a pole in m q , due to 
the exact zero modes. Its residue is 

lim (m q V) 2 (V I (x)P J (0)) u = Tr[T 7 T J ] C H (^)+Tr[T 7 ]Tr[T J ]£ H (z) , (2.65) 
in,/ •() 

where 

K 

C\u\(x-y) = (jTv^ViixW^vjiy))^, (2.66) 

K K 

C\u\(x-y) = ~ (J2 v K*) v i( x )J2 vt M v Ay)) > ( 2 - 67 ) 

i=l j=l 

and the sums are over the set of K — N + + AL zero modes Vj of the Dirac 
operator, T>vi = , V i = 1,...,K. We recall that N + zero modes have 
positive chirality while A_ have negative chirality. The modes are normalized 
so that Jd 4 x vj(x)vi(x) = V. 
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It is important to note that in writing Eqs. ()2.65|) — ()2.67|) it is assumed 
that poles arise only from exact zero modes, i.e. that taking the limit m q — > 
and performing the average over the full space of configurations commute. At 
fixed volume the only potential danger arises from the average distribution 
of eigenvalues near zero; this assumption holds if the density of eigenval- 
ues vanishes for gauge configurations at fixed non-zero topological charge. In 
xPT, as well as in random matrix theory (RMT) [siMs^ . which will be dis- 
cussed in the next Section, the densities behave as p u {\) ~ \( 2 \ v \+ 1+2N f) ; and 
no contribution from the non-zero modes is thus expected in the considered 
observables. 

The zero mode contribution to the pseudo-scalar correlation functions 
(P 1 (x)P I \y)) v at next-to-next-to- leading order (NNLO) in the e-regime for 
a given topological sector v was computed in Ref. [30], both for full and 
quenched QCD. In this computation the chiral expansion in momentum and 
mass was combined with the 1/N C expansion 



P 2 ~ y ~ e 2 • (2-68) 



The corresponding Lagrangian for the quenched xPT to the NNLO reads 

'(2) _ K 

-gxPT 4 



^Spt = ^str [a^tr 1 ] - ^Str [u e u + u- l u^]- 



777 rv 

zm^c&oStr [U U - U^ 1 } + + ^-(^$ ) 2 (2.69) 



It has an additional low energy constant K coupling the singlet scalar field 
$o to the pseudo-scalar meson field U. 

In this particular expansion the mass m was chosen to be a small pa- 
rameter of order e, in contrast to the standard counting where it is treated 
to be of order 0(1). It is also assumed that F n , S and K scale as 

F* = F*N C , S = SiV c , K = ^=. (2.70) 



This implies that 0{F n ) ~ £>(£) ~ 1. Given the fact that 0(m q Y>V) ~ 1 we 
are led to the relation 

m q ~ e 6 . (2.71) 

It is further assumed that 0(a) ~ 0(K) ~ 1 and 0(mfy ~ e 2 . In this 
framework the pseudo-scalar density P 1 corresponding to V 1 in QCD takes 
the form 



P 1 = z-Str 



T 1 (u e U -U^Ug 1 ^ -K$ Str T 1 (UoU+U^Ug 1 ^ . (2.72) 
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The authors of the Ref. worked out the expressions for the first derivative 
of the pseudo-scalar as well as the scalar correlators at non-zero time separa- 
tions in quenched as well as in full x?T. It was pointed out that taking the 
first derivative helps to get rid of the NNLO Gasser-Leutwyler coefficients 
for the expressions in full 

\c? H {t) = Km o (m q V) 2 ±Jd 3 x(P a (x)P a (0)) u , 

N f C' lul (t)+N*C' H (t) = hm Q (m q V) 2 ^Jd 3 x(P (x)P°(0)) u .(2.73) 

C\„\(x) and C\ u \(x) are quantities which correspond to C\ v \{x) and C\ u \(x), 
but calculated in the framework of x?T. If we define the volume average 
C\ u \(t) = J d 3 xC\ u \(x) and C\ v \(t) = J d 3 xC\ v \{x), then at large time t 

C\u\(t) = C H (t) , C\ v \{t) = C\ v \{t) , (2.74) 

which represent the matching of the results in quenched QCD and in 
quenched x?T. 

In the following we redefine a as 

a 2 4^- (2.75) 



2N C 2N C S 

This allows us to simplify the resulting expressions for the zero mode con- 
tributions to the pseudo-scalar correlator in a way that only two low energy 
constants will appear instead of the combination of K, F n , S and a. The 
expressions for qxPT in the volume V = L 3 T read 



2\v\ 



W\K(t) + 



a 



01 



2N C F?VV 
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FZV \ 3 



2N C FZVV 



u\ - 2\u\(u 2 } )h' 2 (r) + \u 



2F%V 
K(r) 



K(r) 



T 2 



'2F3V 



K(r) 



(2.76) 



[2.77) 



where h\{r) and r were defined in Eq. (j2.61|) and /^(t) = 
2^ [r 2 (r — l) 2 — g^]. 0i is a geometrical shape coefficient depending on the 
form of the volume. 0\ as well as hi, h% and h 3 were introduced in Ref. 0] 
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in context of computation of the two-point functions in the a-model at finite 
volume. hs(r) is defined by 



where the prime in the sum indicates that the zero mode term, with ft — , 
is left out. The expectation value (y 2 ) is related to the mass itlq in analogy 
to the Witten-Veneziano formula which will be discussed in Chapter 7. 

2.4 Chiral random matrix theory 

The appearance of the theory of random matrices in the mathematical lit- 
erature goes back to 1928 and was first applied to physics in the context 
of nuclear resonances by Wigner almost 50 years ago. At that time, theo- 
retical approaches such as the shell model had proven to be very successful 
in describing the low-lying excitations of complex nuclei. However, highly 
excited resonances, which can be observed experimentally by neutron scat- 
tering, could not be described by the microscopic theory. The problem is 
generic: for any complex quantum system containing many degrees of free- 
dom with complicated dynamics, it is very hard, if not impossible, to obtain 
exact results for the energy levels far above the ground state of the system. 

Even though the highly excited states cannot yet be predicted individ- 
ually, one can try to find in the experimental data some generic statistical 
features that can be described theoretically. This is where RMT comes in. It 
aims at a statistical description of the energy levels of the system using only 
basic symmetry considerations as its theoretical input. 

In the following discussion we would like to restrict our attention to the 
case of QCD. One can write a Hamilton operator for QCD in matrix form. 
Since QCD is not solved analytically one may follow an approach whereby we 
assume that all interactions which are consistent with the symmetries of the 
system are equally likely. This means that we replace the Hamilton matrix by 
a matrix which is compatible with the global symmetries of QCD and with 
elements which become now uncorrelated and distributed according to some 
distribution of random numbers. To obtain definite results, observables such 
as the level density must then be averaged over the random matrix elements. 
This defines a statistical theory of energy levels, which is known as RMT. 





(2.78) 
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In this Section we will discuss the chiral Random Matrix Theory (xRMT). 
The structure of the xRMT partition function is inspired by the instanton 
gas model. Although it is not necessary, it is instructive to think of the field 
configurations as a superposition of N + instantons and N_ anti-instantons. 
Each isolated instanton or anti-instanton has exactly one fermionic zero mode 
with a definite chirality. In total we have v = N + — iV_ for the index of 
the gauge field configuration. At finite separation of instantons and anti- 
instantons the remaining modes are no longer exact zero modes of the Dirac 
equation and give rise to nonzero overlap matrix elements in the submatrix 
W of the Dirac operator [35j]. 

A model describing the zero mode par t of the QCD partition function in 



the Euclidean volume V is defined by 
Z(m„ S) = £ m t ( I ) e««— » / DW P(W) ft det ( ^ ™ ) 

N>0 N + =0 V 7 J f=l V 1 7 

(2.79) 

N = N + + N_ is the number of zero modes, and the matrix W is a N + x N_ 
matrix. The integration measure DW is given by the Haar measure. The 
binomial factor arises since zero modes of each chirality are treated as inde- 
pendently distributed identical particles. fJ>(N) is the distribution function of 
the total number of zero modes. P(W) denotes the probability distribution 
for field configurations. The structure of the overlap matrix, with off-diagonal 
blocks W and W^ and diagonal blocks equal to the quark masses rrif times 
the identity, is dictated by the chirality of the zero modes, see Eq. (j2.19|) . The 
density N/V of the total number of zero modes is kept fixed as V changes. As 



in the instanton liquid approximation to the QCD partition function 37-39], 
it is considered to be an external parameter. 

We now replace the average over all gauge field configurations with an 
average over Gaussian distributed overlap matrix elements with the distri- 
bution function 

P{W) = exp I j- Tr(WW ] ) J . (2.80) 

(3 is the so-called Dyson index which is defined as the number of independent 
variables per matrix element. To this end the matrix elements of W can be 
either real [f3 = 1, chiral Gaussian Orthogonal Ensemble (%GOE)], complex 
\P = 2, chiral Gaussian Unitary Ensemble (%GUE)], or quaternion real [(3 = 
4, chiral Gaussian Symplectic Ensemble (%GSE)]. In the latter case, the 
eigenvalues of T> are doubly degenerate, and the use of Majorana fermions is 
implemented by replacing the determinant by its square root. 
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We will restrict ourselves to the xRMT partition function in a specific 
topological sector, 

Z^K...,m^)= I M[idet(D + m/ )e-^ Tr ^) , (2.81) 

The parameter X will be later identified with the chiral condensate and N 
with the spacetime volume. For large N we assume that \u\ does not ex- 
ceed a/TV so that, to a good approximation, N + = N/2 for large N. These 
assumptions motivated by the instanton liquid model are, however, not ful- 
filled in the lattice simulations with the Ginsparg- Wilson fermions which are 
the subject of the thesis. In fact one observes zero modes either with posi- 
tive or negative chirality. Therefore an aim of this work is also to prove the 
viability of the predictions by the RMT in our quenched QCD simulations. 

This model inherits the following symmetries of the QCD partition func- 
tion: 

• The Ua(1) symmetry. All eigenvalues of the random matrix Dirac op- 
erator occur in pairs ±iA n , or they are zero. 

• The topological structure of the QCD partition function. The Dirac 
matrix T> has index v = N + — N- . This identifies v as the topological 
charge of the model. 

• The flavor symmetry, which is the same as in QCD. For (3 = 2 it is 

SU(A^) <g> $U(N f ), for (3 = 1 it is SU(2iV», and for (3 = 4 it is SU(jV/) 
(two Majorana flavors count as one Dirac flavor). 

• The chiral symmetry, which is broken spontaneously with a chiral con- 
densate given by 

S = lim lim vrp(0)/A^. (2.83) 

mj ^O N—>oo 

(N is interpreted as the dimensionless volume of spacetime.) The 
symmetry-breaking pattern is SU(A r / ) ® SU(A r / ) -> SU(A r / ), 
SU(2JV» -> Sp(2JV/), and SU(N f ) -> 0(N f ) for = 2, 1, and 4, 
respectively — the same as in QCD |4lN43| . 

The chiral ensembles are, however, only a part of a larger classification 
scheme that also includes ensembles for the description of disordered super- 



conductors 441. In total, 10 different families of random matrix ensembles 
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have been identified. They correspond one-to-one to the Cartan classification 
of symmetric spaces. 

We consider the spectral representation of the xRMT partition function. 
The partition function ()2.8H) is invariant under the transformation 

W -> U ] WV , (2.84) 

where the n x n matrix U and the m x m matrix V are orthogonal matrices 
for {3 = 1, unitary matrices for (3 = 2, and symplectic matrices for /3 = 4. 
Using this property the matrix W can be rewritten as 

W = U ] \V . (2.85) 

Here, A is a diagonal matrix with real positive matrix elements In terms 
of the eigenvalues, the partition function (j2.81j) is given by 

z^k, ...,m Nf) = [nd\ k iA(A 2 )^nAf^ iH+ ^ i e-^^n-/ l ( A i+-/) . 

j k k f 

(2.86) 

where the Vandermonde determinant A is defined by 

A(A 2 ) = H(\l - Xf) . (2.87) 

k<l 

In the followin g we will need the a model representation of the xRMT 
partition function 145^ . One arrives at it by writing the fermion deter- 
minant of Eq. 1)2.81)1 as a Grassmann integral. Then one averages over the 
Gaussian distribution function. The resulting expression is subject to the 
saddle point approximation in the large iV limit. We give the expression for 
the partition function at /3 = 2 



M 

One can see that the partition function for the xRMT in the large N limit 
coincides with the leading order expression of the partition function in 
the e- regime characterized by Eq. 1)2.51)1 [sfil ]. This fact leads to profound 
consequences for the observables calculated in the e-regime. We are going to 
discuss the applications for the so-called microscopic eigenvalue density in 
the next Subsection. 



DU (detuy e^ T ^ MU+M ' u ^ , 

U£V{N f ) 

diag(mi, . . .,m N ) . {21 
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2.4.1 Microscopic spectral properties 

We start off by discussing the Leutwyler-Smilga sum rules j^. The micro- 
scopic spectral density of the Dirac eigenvalues arises naturally from them. 
The Leutwyler-Smilga sum rules are obtained by expanding the partition 
function Z v {mf) in powers of rrif for the expression in the e-regime (|2.51|) . 
and Z^ CD (mf) given in Eq. (j2.10j) . and then equating the coefficients. We 
take the topological sector to be v = 0. Comparing coefficients at 0(m 2 ) we 
obtain 

1 V 2 T/ 2 

A fc >0 k 1 

where the sum is restricted to nonzero positive eigenvalues. The Leutwyler- 
Smilga sum rules can be expressed as an integral over the eigenvalue spectral 
density p(X) and spectral correlation functions, 



p(X)dX 



(2.90) 



V 2 Y? J X 2 AN f 
If we introduce the microscopic variable z, 

z = XVY* , (2.91) 



the integral can be rewritten as 

1 f P(vx)dz 



/ z 2 AN 



f 



(2.92) 



Also according to the Banks-Casher relation the small eigenvalues are 
spaced as 

AA = ^ = ^- < 2 - 93 > 

provided that p(0) /V > 0. These arguments naturally define a scale which is 
suitable for the study of the distribution of individual eigenvalues. For this 
purpose, it is convenient to define the so-called microscopic spectral density 

Ps(z)= hm ^p{4t)= J™ E^- A ^)>Qcd, 

n 

YiVrrif fixed T,Vrrif fixed 

(2.94) 

where the average is taken over the QCD partition function. This limit exists 
if chiral symmetry is broken. Ab initio calculation of this quantity in QCD is 
a challenging task which requires a non-perturbative treatment and could not 
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be solved yet analytically. The microscopic spectral density, however, arises 
as a natural quantity for the xPT in a finite volume since it is defined at fixed 
rrif'EV in the dynamical case. By using the relation between the ;\RT at the 
leading order and the xRMT in the large N limit, we can reduce the problem 



to the calculation of the microscopic spectral density for the %RMT 31-34 . 
The corresponding definition of the microscopic spectral density in xRMT is 
given by 

n 

Y,Nm { fixed EiVm/ fixed 

(2.95) 

where N is the size of the random matrix, (. . . ) x gue denotes averaging with 
respect to the xGUE partition function and A n are the eigenvalues of the 
random matrices. If N corresponds to the lattice volume then E is assumed 
to be taken in dimensionless units. 

The conjecture is that p s (z) is a universal function that only depends on 
the global symmetries of the QCD partition function and can be therefore 
computed by means of the xGUE, 

Ps(z) = P* GVE (z) ■ (2.96) 

Below we illustrate how this statement fits in the framework of xPT and 
X GUE: 

Z? CD (m f ) = I DA^ e~ s W x H^ xUE = / DW e~ Nv ^ ]Jf det(Z> + m f ) 
T\ f det(D + m f ) 

I (a) (c) | 

Zv{M) = I sv{Nf) DU (det U)»x H ...,m Nf ) = 

e Y* Tr {M U+M C7t) = J DW e -^TrWWl ^ de ^ V + m/) 



We show the relation between the partition functions of QCD ()2.9|) . 
%PT (EIHrj and xRMT f2M\ with GUE. 

If we formally substitute the integration over gauge fields in the QCD 
partition function by integration over the Dirac matrices then we come to 
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the formulation in terms of the random matrices W and the partition function 
of the chiral random unitary ensemble (link (b)). We introduce the potential 

v(W) = ^(S YM [A]-\og(J)), (2.97) 

where J is the corresponding Jacobian and N is the size of the random 
matrices, that can be related to the physical volume V. This reduction is, 
however, not well theoretically founded. 

The link (a) corresponds to a reduction of the QCD partition function 
i^QCD in the finite volume to the low energy effective theory for the soft pions, 
which amounts to neglecting the derivative terms in the chiral Lagrangian. 

Once we have a formulation in terms of the unitary random matrices 
Z x tje we can apply the result of Ref. |4f|, where it has been proven that the 
microscopic spectral densities, and all microscopic spectral correlators are 
universal within the given classes of matrix model ensembles. Hence, the link 
(c) leads to the xGUE partition function. 

The link (d) corresponds to the reduction of the yGUE in the large N 
limit to the partition function for the soft pions [36j]. 

All this pattern can be reformulated for quenched QCD. Since RMT is 
well defined for any number of flavors, one needs only to substitute the 
by its quenched formulation. Also the formalism of the qxPT is known. 

There is now mounting evidence that the described scenario is correct. 
The microscopic spectral densities, derived from large N random matrix en- 
sembles with Gaussian weights |47H5l| . have been shown to consistently re- 
produce the exact spectral sum rules of Leutwyler and Smilga. However there 
is no similar analytical check for the probability distributions of the individ- 
ual eigenvalues. 

In particular we want to highlight the analytical results for the micro- 
scopic spectral density and the probability distributions of the individual 
eigenvalues. The latter were calculated in the framework of the xRMT with 
the Gaussian unitary chiral ensemble. 

The microscopic spectral density was calculated for the xGUE, as well 
as for the e-regime of quenched x?T. (Also some other s pec ial cases were 
computed in \PT which agree well with %RMT results |48M5ol The 
microscopic spectral density in the topological sector \v\ for the xGUE reads 

P { s\ z ) = \ {j 2 N f + \v\{ z ) ~ JN f + \u\+l{z)JN f + \v\-l{z)} + \v\6(z) , (2.98) 

where the J's are Bessel functions. These expressions can be reassembled in 
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the microscopic spectral density p s (z), 

oo 

p s (z) = P^s\ z ) ? v — topological charge . (2.99) 

v=— oo 

These functions pf\z) can be expressed in series of distributions of the k th 
individual lowest eigenvalue for each topological sector. Refs. jHHIH present 
the general formula for the probability distributions of the individual eigen- 
values Pk(z) calculated for the xGUE, 

^*) = £^*) + M*(*) > (2-ioo) 

fc>l 

where the sum runs over the non-zero modes. The distribution of the lowest 
eigenvalue of the Dirac operator for an arbitrary topological sector and Nf = 

A V \z) = |e- z2 /Met[/ 2+l _^)] M=1 ,..., H . (2.101) 

The predictions for the lowest eigenvalue at \v\ = 0, 1 and 2 are depicted 
in Figure 12.11 We see in particular that the density peak moves to larger 
values of z as \u\ increases. 




Figure 2.1: The distributions for the lowest eigenvalue at \v\ = 0, 1 and 2, as 
predicted by xRMT in Eq. (jruTT^ . 
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2.4.2 Bulk eigenvalues. Unfolding of the spectra 

The spectra of many complex quantum systems exhibit a universal behavior. 
To see this, however, one needs to unfold the spectrum. The unfolding can 
be done if the scale of the variations of the average spectral density is much 
larger than the scale of the spectral fluctuations. Then we can separate the 
spectral density in an average density p(A) and a fluctuating part, 



The details of the separation are specific for the observable of interest. The 
only common feature is that the resulting unfolded spectrum {x n } has aver- 
age spectral density equal to 1 and is given by 



where {E n } represents the original spectrum. There are several ways to ob- 
tain the average spectral density. In some cases it can be conjectured analyt- 
ically by some polynomial approximations. In most cases, however, we have 
to average over many level spacings. There are two different procedures to 
achieve this. The first is the spectral averaging where one takes averages over 
many level spacing around the energy of interest. The second is the ensem- 
ble average where one performs averages over many different level spacings 
drawn from all configurations of the ensemble. If we calculate the spectral 
correlation functions from the unfolded spectrum then these two procedures 
do not necessarily yield the same behavior at long distances. The equivalence 
of the two procedures with respect to the correlation functions is known as 
"spectral ergodicity" . 

The distribution of the unfolded eigenvalues takes the usual form 



In this context one usually considers the unfolded level spacing A n as a 
tool to probe the short range correlations in the spectrum 



p(A)=p(A)+ Pfl (A) . 



(2.102) 




(2.103) 




(2.104) 



n 



X n X. 



n—1 



(2.105) 



Its distribution is given by 




(2.106) 



n 
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This function is equal to the probability density for two neighboring levels n 
and ri + 1 having the spacing s. 

For the uncorrected case it is expected that the level spacings are dis- 
tributed according to the Poisson distribution p(s) = exp(— s). A good ansatz 
for the correlated case is provided by the Wigner distributions, which are the 
exact "spacing distribution" for a 2 x 2 matrix model 



_7T 2 
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s e f 



262144 „4 



orthogonal ensemble 



unitary ensemble 



(2.107) 



< 729tt 3 



s e symplectic ensemble . 



2.5 Lattice QCD 

2.5.1 The Wilson gauge action and the Wilson 
fermions 

QCD is a renormalizable quantum field theory that must be regularized. A 
way to achieve this is to put the system on a cubic space-time lattice with 
some lattice spacing a jH^. This introduces a momentum cut-off oin/a. Now 
one has to introduce gauge and quark degrees of freedom on the lattice. 
Maintaining manifest gauge invariance is essential when gauge theories are 
regularized on the lattice. In order to achieve this one needs a convenient 
way to represent the gauge fields on the lattice. Wegner and Wilson, as 
well as Smit, independently introduced the concept of a parallel transporter 
U X)IX G SU(N C ) connecting neighboring lattice points x and x + afi. Here fi is 
a unit vector pointing in the /i direction. The parallel transporter is related 
to an underlying continuum gauge field A^(x) = igA a ^(x)T a by 

U x> „ = Pexp [ ds A„(x + fls) , (2.108) 
Jo 

where V denotes path-ordering for non-Abelian gauge theory. Under a gauge 
transformation Q x the parallel transporter transforms as 

= a-r,,, <>,.,, . (2.109) 

One tries to preserve as many symmetries as possible at the finite lattice 
spacing. The hope is that the rest will be recovered in the continuum limit. 
A gauge action that respects the gauge symmetry was first put forward by 
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(p) 

K. G. Wilson. It is a sum over the so-called plaquette variables U X J V , which 
are the products of link variables around a closed square path, 



fiuJ i 



Sp(Ui%) = 1-1-ReTrU^, 

uiju = ^x,nU x +^ jV U x ^. Vtft U^. >v . (2.110) 

U x +^ v represents a link attached to a site x + afi and pointing in the direction 
v and Ul v represents a link pointing in the opposite direction to v and 
attached to the site x + av. [3 is related to the bare coupling constant as 
P = ^r- The Wilson gauge action reduces to the continuum Yang-Mills 
action (|2.5J ) in the naive continuum limit a —>■ 0. This is demonstrated e.g. 
in Refs. 0,1^. 

To fully define the path integral we must also consider the measure. The 
lattice functional integral is obtained as an integral over all configurations of 
parallel transporters U Xj/Ji , i.e. 



n 



X,fl 



[ dU x>ti ex P (-S YMW [U]) . (2.111) 

JSU(N C ) 



One integrates independently over all link variables using locally the Haar 
measure dU^ x for each parallel transporter. 

For compact groups like SU (N c ) the integration is limited to a finite do- 
main. This makes it unnecessary to fix the gauge in lattice QCD because the 
functional integral is finite even without gauge fixing. This is an important 
advantage of the formulation using parallel transporters. 

The Yang-Mills functional integral contains a single parameter — the 
bare gauge coupling g. The continuum limit is taken by searching for values 
of g for which the correlation length of the lattice theory diverges in lattice 
units. In the language of statistical mechanics one looks for a second order 
phase transition. Due to asymptotic freedom, in lattice QCD one expects a 
second order phase transition at g — *■ 0. 

The physical observables must be formulated in a gauge invariant way. In 
fact, the path integral average over gauge non- invariant observables vanish. 
This statement is a consequence of Elitzur's theorem |Hf| which in particular 
implies 

(U x>fl ) = (2.112) 

for each link variable U Xttl . 
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Now we would like to discuss how the fermion fields are put on the lattice. 
The quark lattice fields \l/ x are Grassmann variables attached to the sites 
x. The Dirac operator contains a continuum derivative which has to be regu- 
larized on the lattice. It turned out that it is a highly non-trivial problem to 
regularize it in a way that preserves the chiral symmetry. The latter becomes 
of great importance when addressing the physics of light quarks. 

A simple way to regularize the Dirac operator is the so-called naive dis- 
cretization of fermions on a lattice. One substitutes the partial derivatives 
with finite differences, 

*] = a 4 ^x{D) xy ^ y = a 4 ^ — (^7^+£-^+£7^)+ a4 Yl m ^^* 

x,y x,fi x 

(2.113) 

This action has the so-called doubling problem which is also present in the 
interacting case, where the action reads, 

X,fl X 

The free quark propagator has poles at momenta p M = {0, ir /a}, which gives 
rise to a total of 2 4 species in 4 dimensions. 

It was Wilson who proposed to eliminate the unwanted doubler fermions 
by breaking chiral symmetry explicitly [5?}. The unwanted 2 4 — 1 species can 
be sent to the cut-off scale by adding to the Dirac operator an additional 
term, the Wilson term, which removes the doublers in the continuum limit, 



x,y 

X,fi x 

+ a 4 £ ^(2**** - ^ x U x ,^ x+(l - * x +pPl^ x ) ■ (2.114) 



x,ii 



The doubling problem is expressed by the Nielsen and Ninomiya No-Go the- 
orem 58^ 60| which states that it is not possible to simultaneously solve the 
doubling problem and have exact chiral symmetry on the lattice with a local, 
discretely translation invariant and real Euclidean action for lattice fermions. 
The locality of a lattice action is defined by requiring the couplings to decay 
exponentially with the distance. This provides a safe continuum limit for a 
lattice Dirac operator. 
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Many lattice Dirac operators, such as the Wilson Dirac operator, satisfy 
the following condition, which is called 75-Hermiticity, 

& = 75^75 • (2.115) 

If the lattice is finite, one can simulate the theory on a computer. This can be 
done efficiently only in Euclidean space, where the gluonic weight function 
is exp(— SVm) with SVm real and positive. The full weight function of QCD 
also contains the fermion determinant, which can be expressed in terms of 
the gauge fields. Observables are computed by generating gauge field configu- 
rations at some value of (3 in a Monte Carlo update procedure and averaging 
an observable over many configurations. 

Since the inclusion of the fermion determinant is very time consuming, it 
is a well-known simplification to use only the gluonic part of the weight func- 
tion in the Monte Carlo updates, i.e. to perform the quenched simulations, 
which correspond in some sense to the limit Nf = 0. To make contact with 
continuum physics, the results of lattice simulations must be extrapolated to 
the continuum limit, i.e to the infinite lattice volume and zero lattice spacing 
keeping the physical volume fixed. 

2.5.2 Staggered fermions 

For the sake of the following discussion we would like to introduce in this Sub- 
section the staggered fermions. Staggered fermions are obtained from naive 
doubled lattice fermions (|2. 113(1 by the so-called spin diagonalization 
It consists in performing a local change of fermionic variables 

* t -»4*,, (2.116) 

47^+aA = A m(*) e U(lf d . (2.117) 

A x is a d x d unitary matrix which performs the spin diagonalization and 
A M (x) is a diagonal unitary matrix belonging to the direct product of U(l) 
groups denoted by £7(1). In four space-time dimensions the size of the Dirac 
matrices is 4 x 4. By spin diagonalization one can reduce the fermion multipli- 
cation factor 2 d to 2 d /4, where d is the spacetime dimension. Hence, for d = 4, 
staggered fermions represent 2 4 /4 = 4 pseudo-flavors of mass-degenerate 
fermions. The new fermionic variables \ x become a 3 (colors) x 1 (spin) com- 
ponent object. The corresponding lattice action for free staggered fermions 
takes the form 

S[x, x\ = ^ ^^{XxVx^Xx+fi - Xx+^x^Xx) + a 4 Yl m XxXx , (2.118) 

X,fl X 
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where 

Th,l = h V*,2 = (-l)^ , tfc,3 = (-1) ( - X1+X ^ a , VxA = (_l)(- 1 +-2+X3)/a 

(2.119) 

The 4 spin degrees of freedom for the 4 pseudo-flavors are made from the 16 
components of the x field on a 2 4 hypercube 

For m = the interacting staggered fermion action has an exact remnant 
[/(l) e ®t/(l)o chiral symmetry 

X x exp(— iipo) , for {x\ + x 2 + x 3 + x±)/a even , 
Xexp(— i<f e ) , for (x\ + x 2 + x 3 + x^j/a odd , 

(2.120) 

which is a subgroup of the SU(4)l <S> SU(4)r ® ?7(1)b chiral symmetry of 
the corresponding continuum theory. These symmetries are expected to be 
recovered in the continuum limit. 

Due to the remnant chiral symmetry the staggered fermions are protected 
from some problems which are inherent for the Wilson fermions. These are 
the additive mass renormalization and the problem of the exceptional config- 
urations in the quenched simulations, i.e. the occurrence of the configurations 
with accidentally small eigenvalues of the lattice Dirac operator. The stag- 
gered fermions are also order a 2 improved i.e. the naive discretization errors 
are expected to behave like 0(a 2 ) which is opposed to the Wilson fermions 
where the errors are of order 0(a). 

The first down-side of staggered quarks is the fact that we have to deal 
with exactly four pseudo-flavors. Secondly, the discretization errors induce 
flavor-changing interactions and so are rather dangerous. To overcome the 
first problem it is tempting to consider a fourth square root of the staggered 
Dirac operator to obtain only one pseudo-flavor theory. This operation, how- 
ever, is dangerous since the locality of the resulting action is questionable • 

2.5.3 The Ginsparg— Wilson relation and the Neu- 
berger overlap operator 

It were Ginsparg and Wilson who first put forward a relation on the lattice 
Dirac operator which helps to circumvent the Nielsen and Ninomiya No-Go 
theorem. One gives up the chiral symmetry in its conventional formulation, 

{A 75} = 0, (2.121) 



66 . 



X x = exp(z^ e )xx , X x = 
Xx = exp(^ )Xz , % = 



38 



Chapter 2 Theoretical background 



and requires for the lattice Dirac operator to satisfy the remnant chiral sym- 
metry condition of Ginsparg and Wilson [68^ . 

D l5 + l5 D = -D l5 RD . (2.122) 
H 

where R x , y is a local operator for the x, y indices and its structure in the 
spinor space is dictated by {75, R} = 0. The Ginsparg- Wilson relation states 
that the continuum condition for the right-hand-side to vanish is relaxed to 
a term of 0(a), where a is the lattice spacing. \i is a mass parameter of order 
(9(1), see below. The Ginsparg- Wilson relation (GWR) implies that 

ZTSs + = (a/tih 5 R (2.123) 

is local. It is a formulation of the chiral symmetry after performing the Wilson 
renormalization group transformation from the continuum theory to the lat- 
tice. Hence the Ginsparg- Wilson fermions necessarily have exact zero modes, 
which occur with positive or negative chirality. From now on we will consider 
Rxy = S xy . We follow the idea of M. Liischer and write down a modified 
infinitesimal chiral transformation on the lattice 

= + tftf = M + ie a T a j 5 (l - ^D)j * > 
^' = ^ + 5^ = ^ (l+ie a T a (l-^-D)^ . (2.124) 



where e a are parameters of the chiral rotation |6£ 

Through D Liischer's lattice version of a chiral transformation depends on 
the gluon field. Still, in the continuum limit a — > it reduces to the standard 
chiral symmetry of the continuum theory. Any lattice fermion action of the 
form Sp = a 4 ^Dty is invariant under transformation (|2.124j) to order 0(e), 
provided that D obeys the Ginsparg- Wilson relation Eq. ()2.122|) . Indeed, 



a^Die a T a lh (l - —D) + a A ^ie a T a (l - —D)y 5 DS> 



SS F = a\6VDy + yD6y) 
i^Die a T a l5 (l - £ 

= a 4 me a T a {^ 5 , D}^ - a 4 ^ie a T a -D^Dm = . (2.125) 

Under this transformation the fermion measure D^Dty is not U(1)a in- 
variant. This is analogous to continuum QCD where this property of the 
fermion measure in the partition function leads to the anomaly for the axial 
current. In this way the axial anomaly is also correctly reproduced by the 
Ginsparg- Wilson fermions 69T 7l|. 
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We now define the topological charge v simply by the index 



N:~N_ 



(2.126) 



where N + (iV_) is the number of positive (negative) chiral zero modes. Here 
one adopts the continuum index theorem and uses it to define the topological 
charge of the lattice gauge configurations Since the Wilson operator 

D w obeys 75-Hermiticity f)2.1 15j) this allows for a simple solution of the GWR 
by inserting D\y into the overlap formula 



a 



1 + A/VAA 



A = aD 



(2.127) 



D = D 



w 



(2.128) 



D ov is a solution to condition ()2.122|) and it is 75-Hermitian as well, fi rep- 
resents a negative mass of the Wilson fermion, which can be chosen in some 
interval as long as the gauge fields are smooth. This particular solution with 
the kernel D = Dw is often referred to as the Neuberger overlap operator. 
In this thesis we will denote it also as the overlap Wilson Dirac operator. 

If we assume that the operator D already satisfies the GWR and also is 
75-Hermitian, then 



A ] A = (at) ] - n)(aD - y) = a 2 D ] D - u.a(D ] + D) + a 2 



r/7/7-, ( ^Dnf S b - {A 75} 1 + /' 



2 



(2.129) 



This implies for the overlap operator D c 



a 



1 + 



aD — n 



D 



(2.130) 



i.e. an operator satisfying the GWR is reproduced by the overlap formula [7 

Note that the operator A/y MA is unitary, hence the spectrum of D ov is 
located on a circle in the complex plane through zero, with center and radius 
fi/a, which we denote as the GW circle. 



2.5.4 The hypercube Dirac operator 

There is a natural way to proceed from the continuum formulation of the 
action to a lattice form. In this approach one defines the lattice fermion 
fields as block averages of continuum fields integrated over hypercubes. The 
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lattice action of these lattice fermion fields is related to the continuum one 
by the Wilson renormalization group transformation. The resulting lattice 
theory, which is called the perfect lattice fermions, is equivalent to the un- 
derlying continuum theory with respect to the physical observables, i.e. it 
is completely free of lattice artifacts. Similarly one can also start from a 
given lattice formulation on a fine lattice and — by using the renormaliza- 
tion group transformation — arrive at another lattice formulation at coarser 
lattice which is physically equivalent. Iterating this procedure infinitely many 
times we arrive at the fixed point action (here we assume that this limit ex- 
ists). The fixed point action of the quantum renormalization group trajectory 
represents perfect lattice fermions since it is insensitive to a change of the 
lattice spacing. 

To simplify the notation we will set the lattice spacing to a = 1. The 
hypercube Dirac operator is a free fixed point action with its couplings being 
truncated to the unit hypercube. In this Subsection we consider the con- 
struction of the hypercube Dirac operator. The discussion will closely follow 
Ref. jlfj. We will start from the case of free fermions and then explain how 
gauge fields can be introduced. 

Let us divide the (infinite) lattice into disjoint hypercubic blocks of n d 
sites each and introduce new variables living on the centers of these blocks 
(block factor n RGT). Then the RGT relates 



where \1/ represents the fermions on the fine lattice and the represents 
fermions on the coarser lattice. The points x G are the sites of the original 
fine lattice and x' are those of the new lattice with spacing n. x G x 1 means 
that the site x belongs to the block with center x'. 

Now the original action S[$, transforms into a new action S'[§?', \&'] 
on the coarse lattice. The latter is determined by the functional integral 



The kernel K[$', ty', has to be chosen such that the partition function 
and all expectation values remain invariant under the RGT. At the end, one 
usually rescales the lattice spacing back to 1. In any case, the correlation 
length in lattice units gets divided by n. 

For the kernel functional there are many possible choices j??], [7^ . We will 
use the Gaussian type kernel 
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x ex p { - ^ E - E *«] - ^1)7^ E *$ 133) 

This type of transformation with non-vanishing a > is not chirally invari- 
ant, which will lead to the formulation with broken chiral symmetry in its 
standard form given by Eq. (j2.121|) . 

Assume that we are on a "critical surface", where the correlation length 
is infinite. After an infinite number of RGT iterations we obtain a finite fixed 
point action (FPA) S 1 *^, The critical surface requires a fermion mass 
m — 0, but one can generalize the consideration to a finite mass. 

Assume that we want to perform a number N of RGT iterations. If we 
start from a small mass mj (nN), then the final mass will be m. In the limit 
nN — > oo (i.e. we start from an infinitesimal mass) we obtain a perfect 
action at finite mass. In this context, "perfect" means that dimensionless 
quantities do not depend on the lattice spacing, hence they are identical to 
the continuum values. 

For the above transformation ()2.133|) . this perfect action can be computed 
analytically in momentum space |Z2(. The computation simplifies if we let 
n — > oo, so that iV = 1 is sufficient. Hence we start from the continuum 
action now, and the perfect action takes the form 

= f d d p ^(-p)A*(p)-^(p) , 

2sin(p M /2) 



= ^ ! (2134) 
;.i ft- 

where A* is the free perfect propagator. The same perfect action is obtained 
starting from a variety of lattice actions in particular from the Wilson 
fermion action. 

In coordinate space we write this action as 

S*[%*] = E^[P»7„ + A(r)]* s+r . (2.135) 

x,r 

At the cost of broken chiral symmetry we obtain a theory with its cou- 
plings in and A decaying exponentially as |r| increases. An exception is 
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the case d = 1, where they are confined to one lattice spacing for the special 
choice 

e m — m — 1 , 
a = = . 2.136 

m z 

It is also in agreement with the Nielsen-Nino miya theorem that, in order 
to obtain a local perfect action, one must break chiral symmetry explicitly. 
Although the chiral symmetry is not manifest in the perfect action, all chiral 
properties are still correctly reproduced by it. This is due to famous Ginsparg- 
Wilson relation, introduced before in (|2.122|) . to which the perfect action 
obeys, 

{A*(p), 75 } = 2a 75 . (2.137) 
Here we can identify a with the mass parameter //, namely, 

a = -L . (2.138) 

It turns out that for the choice of a in Eq. ()2.136|) the locality is also 
excellent in higher dimensions, i.e., the exponential decay of the couplings is 
very fast. This is important, because for practical purposes the couplings have 
to be truncated to a short range, and the truncation should not distort the 
perfect properties too much. Ref. [^H proposed a nice truncation scheme that 
uses periodic boundary conditions over 3 lattice spacings and thus confines 
the couplings to a unit hypercube. It was pointed out that the spectral and 
thermodynamic pro per ties of the HF are still drastically improved compared 
to Wilson fermions jsii li^ . 

It is far more difficult to construct an approximately perfect action for 
a complicated interacting theory like QCD. However, as it was proposed 



in Ref. 76, 81, 83] one can just use a simple ansatz for HF together with 
the standard gauge link variables. Apart from nearest neighbors, one also has 
couplings over 2, 3 and 4-space diagonals in the unit hypercube. One connects 
all these coupled sites by all possible shortest lattice paths, by multiplying the 
compact gauge fields on the path links. This procedure was called "minimal 
gauging". Note that one can connect two sites x and y lying on <i-space 
diagonals via d\ such shortest lattice paths. One averages over all of them to 
construct the hyper-link, see Figure I2~2l which we adopted from Ref. j76j]. The 
hyper-link Up (x) between site x and x + fx is identified with U Xjjjl , and we 
denote the hyper-link in plane, cube and hyper-cube as U^ u (x), U^_(_ v+p {x), 

and Uj^. v+p+a (x) , respectively. Then we can write the corresponding fermion 
matrix in terms of the hyper-links which are constructed recursively starting 
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from the gauge links U f 




(2.139) 



Figure 2.2: 1-space, 2-space and 3-space hyper-links. 

It is convenient to introduce prefactors which are functions of the HF 
hopping parameters Ki and Aj, i = 1, . . . , 4, and sums of 7-matrices 

r± M = Ai + /«i(±7 M ) 

= A 3 + k 3 (±7 m ±7^ ±7 P ) 
ri^ztvipia = A 4 + k 4 (±7m ±7i/ ± 7p ± 7a) • (2.140) 

Note that the Aj in Eq. (|2.14()jl differ from A(r) in Eq. (|2. 135)1 by a normal- 
ization factor The ft, arise from p M (r) by the same normalization. The 
couplings are given in Ref. 0] and in Chapter 4. 

The hypercube operator is organized in sums which run over four different 
directions for two 1-space links, six directions for four 2-space links, four 
directions for eight 3-space links, and one direction for the sixteen 4-space 
links. Altogether 80 hyper-links plus the term 5 xy contribute. With each path 
of the free HF, a 7-matrix is associated. The 7-matrices are chosen so that 
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they add up to produce a prefactor T, see Eq. (|2.140jl . which is associated 
with a given hyper-link. 

The 1-space links llj^\x) are identical to the Hermitian-conjugate links 
in negative direction, Ul + ^ , and this feature also holds for the hyper-links, 
e.g. 

^S M2+ ... +w (x) = U^^—J? + h + ^ + • • • + Ad) • (2.141) 

Therefore, only one half of the 3 d — 1 hyper-links has to be computed and 
stored in the implementation of the HF. 

As in the case of Wilson fermions, the HF matrix exhibits the "75- 
Hermiticity" 

1 5 D HF75 = D^ HF , (2.142) 

i.e. Djjp is non-Hermitian but its eigenvalues come in complex-conjugate 
pairs. The corresponding HF matrix is defined by 



45 



Chapter 2 Theoretical background 



D HF (x,y) 
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2.5.5 The overlap hypercube Dirac operator 

In this Section we would like to present motivations for a non-standard choice 
of the kernel in the overlap formula. 

The generalization of the overlap fermions to a whole class of solutions 
of the GWR and the motivation for considering alternative kernels D were 
given in Ref . jz^| . In particular it was pointed out that the Ginsparg- Wilson 
fermions reproduce themselves if the kernel already obeys GWR as we showed 
in Section 2.5.3. To put it the other way: if the kernel does not obey the GWR 
then it gets corrected by the overlap formula 75|. This additional freedom of 



having different kernels can be used in order to improve such properties of the 
overlap operator as locality and get rid of the lattice artifacts which reside, 
for example, in approximate rotational symmetry. Also the dependence of 
the observables on different lattice spacings, that is referred to as the scaling 
behavior, is expected to be improvable. To this end we start from the kernel 
that already has good properties in rotational invariance, locality, scaling 
behavior and approximately satisfies the GWR. If we insert it into the overlap 
formula then we can expect that due to its already good chiral properties the 
chiral corrections to the operator will be small and therefore the properties 
of locality, approximate rotational symmetry and scaling behavior will be 
essentially inherited by the resulting overlap operator |7a |. A good candidate 
for that would be the approximate perfect fermions introduced in the previous 
Subsection. 

When constructing the overlap operator the parameter /i must be chosen 
in some interval to ensure that the doublers are projected to the right arc of 
the GW spectrum. This interval is fi G [0, 2] for the free Wilson fermions and 
shrinks from both sides in the interacting case. This procedure can be realized 
once the gauge configurations are smooth enough (or (3 is sufficiently large) to 
allow for the spectrum of the kernel develop such a window. It is conceivable 
that using the hypercube operator as a kernel in the overlap formula one can 
admit lower values of f3 than it is the case for the Wilson kernel. This is based 
on the fact that the hypercube operator approximately obeys GW relation 
and this should lead to a better behavior in its spectrum. Having the doublers 
projected to the cut-off scale one can optimize the overlap operator further 
by requiring the best locality properties or the least condition number for 
the A* A operator. This will specify /i and also impose additional constraints 
on the couplings in the hypercube operator. 
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3.1 Quenched simulations of the gauge fields 

The QCD partition function in the formalism of the path integral is defined 
by an integral with an infinite number of dimensions which cannot be ana- 
lytically solved so far. However its regularization on a finite lattice is reduced 
to an integral over a finite number of dimensions and can be used to address 
non-perturbative physics. To this end one can use the Monte Carlo (MC) 
simulations to evaluate the lattice regularized partition function. 

In lattice MC simulations of quenched QCD one usually generates equilib- 
rium gauge configurations of SU (3) fields U with a probability distribution 
proportional to VK[C/] = e~ SvM ^ at a certain value of the strong gauge 
coupling. This procedure is referred to as importance sampling. It can be 
simulated by generating a Markov chain of configurations. One parameter- 
izes the Markov chain with r and introduces the transition probability matrix 
P{\U T+ i] <— [U T ]) connecting the field configuration at simulation time r and 
t + 1. A condition to arrive at the equilibrium distribution in the limit of 
large r is detailed balance 

P( [U T+1 ] <- [U T ] ) W[U T ] = P( [U T ] <- [U T+ i] ) W[U T+1 ] . (3.1) 

For quenched simulations it is numerically favorable to perform a link by 
link update of the given configuration instead of the update of it as a whole. 
The process when all links of the field configuration have been once updated 
is referred to as a configuration update or one sweep. In quenched £77(3) 
gauge simulations two transition probability distributions are mostly used. 
These are the heat bath js^l and Metropolis js^J transition probability dis- 
tributions. The heat bath was conventionally designed for the generation of 
SU (2) pure gauge configurations. To adopt it for quenched simulations with 
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the 577(3) gauge group one needs to use the Cabbibo-Marinari update 
whereby one splits an SU (3) link variable in three 577(2) link variables and 
performs the heat bath update over each of them. Then these three SU(2) 
link variables are reassembled again in the updated S77(3) link variable. In 
Markov processes one has to consider the autocorrelation of the configura- 
tions. Indeed since every next configuration is produced from the previous 
one, it contains information about it. This spoils the sampling of observables 
which in the spirit of MC simulations are supposed to be computed on an un- 
correlated set of configurations picked out with the probability distribution 
W[U] = e~ SyM ^. Since the quenched updates are computationally quite fast 
we can arrange for many of them in-between to make sure that the result- 
ing configurations become uncorrelated. To cope with the autocorrelations 
and thereby to speed up the updating process it was proposed also 
to choose the trial elements for the update in a particular way. This is the 
so-called overrelaxation update. It is achieved by choosing the new link U' x 
as far as possible from U x>(i . The overrelaxation update for SU(N) gauge 
fields relies on the particular form of the lattice gauge action. The part of 
the gauge action containing one link variable U x ^ must have the form 

AS[U X J = const ■ ReTr^V) , (3.2) 

where the matrix V does not necessarily belong to SU(N). This is the case 
in particular for the Wilson gauge action ()2.11()j) that we are using. 

3.2 A numerical treatment of the overlap 
Dirac operator 

In this Subsection we discuss numerical details of simulations of the overlap 
Dirac operator. 

The formula of the overlap Dirac operator contains the square root of a 
Hermitian operator. It is instructive to rewrite it as 



1 + 75" 



Dov = A 4 



Q = c 75 (/J-/i) , D = D W . (3.3) 

Q is a Hermitian operator and constant c was chosen to ensure that the 
spectrum of Q 2 is bounded by 1. For our treatment of the square root we 
used the Chebyshev polynomials T n (x) to approximate the function 1/y/x 



49 



Chapter 3 Lattice simulations 



in the interval (e,l), where e is the inverse condition number of the operator 
Q 2 . This method was proposed in Ref. [oil ] 

The Chebyshev polynomials Tj are the orthonormal polynomials with 
respect to the inner product 

(Ti, Tj) = -^== Tl (t)T,(t) dt . (3.4) 

We thus have 

(T h Tj) = 8q • (3.5) 

Every function / for which (/, /) exists can be expanded into its Chebyshev 
series 

oo 

f( x ) = /jCiTjjx) where q = (/, Tj) . (3.6) 

i=0 

Truncating the series at the iV-th summand gives the polynomial approx- 
imation of degree N — 1 

N-l 

Pn-i(x) = ^qTj(x) . (3.7) 

i=0 

The Chebyshev polynomials T n (x) have exactly n zeros on the interval 
[—1,1], and they are located at the points 

/tt(£;-0.5)\ , . , 

Xk = cos where k — 1, 2, . . . , n . (3.8) 

V n J 

In the same interval there are n + 1 extrema (maxima and minima), 
located at 

Xk = cos ^ — J where k — 0, 1, . . . , n . (3.9) 

At all of the maxima we have T n (x) = 1, while at all of the minima T n (x) = 
— 1. This property makes the Chebyshev polynomials so useful in polynomial 
approximation of functions. It provides the basis for deducing the bounds of 
the maximal error in the approximation. 

The coefficients Cj of the series ()3.7j) are given by 

2 N 

c i = jv ^2fM T sM , wnere 3 = 0> ...,JV- 1 , (3.10) 

fc=i 
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where Xk are the N zeros of T^{x). Then for an arbitrary function / on the 
interval [—1,1] the approximation formula 



f(x) ~ P N -!(x) 



'AT-1 



c fc T fc (x) 



fc=0 



2^o (3.11) 



is exact at x\, . . . , a; at. 

For any point x G [—1,1] this gives the Chebyshev approximation of 
degree N — l for function f(x). The Chebyshev approximation is not the best 
approximation with respect to all polynomials of the same degree, i.e. it does 
not have the smallest maximum deviation from the true general function. The 
latter is true for the so-called minimax polynomials. However the Chebyshev 
approximation is almost as good as the minimax polynomials. The advantage 
of the Chebyshev polynomials is that the formulae ()3.10|) and (I3.11|) can be 
arranged in one recurrence process, the Clenshaw recurrence formula [92| . 
This property makes the Chebyshev polynomials attractive for numerical 
simulations since one can avoid the calculation and use of the Cj coefficients 
which can take large magnitudes and bring in numerical instabilities. 

In our simulations we made the following approximation 

-j= = P Nte (t), where t G [e, 1] , (3.12) 

and e = A min /A max , where A min and A max are the minimal and maximal eigen- 
values of Q 2 , respectively. 

Inserting the approximation to the square root we obtain for the overlap 
Dirac operator 

D ov ~fM[l + l5 QP N>e (Q 2 )] . (3.13) 

The degree N of the polynomial approximation was chosen in a way to 
obtain a specific precision e over i a p for the square root approximation 

\\x-q 2 p n aq 2 ) 2 x\\ 

S eoveriap , (A 14 J 



4||X|| 

where X is a random spinor vector and e over i ap was usually chosen between 
10~ 12 and 10~ 15 . We checked that this condition is insensitive to different 
choices of the random vector X. It also implies a similar precision for the 
GWR. 

The computational cost of the approximation is roughly proportional to 
the polynomial degree. The latter is proportional to the square root of the 
condition number. In order to decrease the condition number of the Q 2 op- 
erator, A max /A min , and thus to decrease the required polynomial degree, we 
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projected out the lowest eigenmodes, 

rn 

Qlroi = Q 2 -J2 X l\ k )( k \> where Q 2 \k) = \\\k). (3.15) 



fc=i 

Then the overlap formula takes the form 



D ov = fj, 



1 + 75Qproj-P/V,e(<2p ro j) + 75 ^ \ k ) ( k \ si S n ( A * 



k=l 



(3.16) 



and the condition number has to be evaluated with respect to the operator 
<5p ro j. To render the chiral symmetry reliable and assure the convergence of 
the polynomial approximation, the eigenvectors of the operator Q 2 should 
be computed to a very good accuracy. We take it again to be between 10~ 12 
and 10" 15 . 

For (3 = 5.85 and (3 = 6 on lattices of V — 8 4 , 10 4 we projected out 
10 eigenvectors. It turned out that projecting out more of them does not 
help much any more. This is due to the fact that above a certain energy 
the eigenvalues are densely packed and further increasing the number of 
projected out eigenvectors does not influence much the condition number. 
In Figure 13.11 we show the condition number for the lattice volumes V = 
8 4 and 10 4 at [3 = 5.85 and /i = 1.6. Here we projected out 5, 10 and 
20 eigenvectors, respectively. The corresponding degree of the polynomial 
approximation varies for different lattice volumes and j3. For instance on the 
10 4 lattice at f3 = 5.85, for 44 projected out eigenvectors it was about 180 at 

e i n— 15 

toverlap J- u 

In addition one has to assure that one deals with a local formulation of 
the overlap fermions. The locality of the overlap Wilson fermions is related 
to the spectrum of its Q 2 operator. Let the spectrum be bounded by u from 
below and by 1 from above 

u < X(Q 2 ) < 1 . (3.17) 

If u comes too close to the origin and the spectrum develops also very dense 
eigenvalues then this will lead to strong fluctuations in the inverse square 
root and locality will be lost. Note, however, that the presence of a few 
near-zero eigenvalues does not spoil the locality. The issue of locality was 
addressed in Ref. [9l| where also a numerical study at /3 = 6, 6.2, 6.4 was 
presented. Analytically the locality properties of the overlap operator can 
be established for smooth configurations. Once each plaquette variable 
defined in Eq. 1)2.110)1 obeys an inequality of the form 

1 1 - ^Re Tr U {p) \<e for all plaquettes (3.18) 
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Figure 3.1: The condition number of the Q 2 operator in V = 8 4 and 10 4 
volumes at (3 = 5.85 and \x = 1.6 after m = 5, 10 and 20 eigenvectors being 
projected out. 
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then the operator Q 2 is uniformly bounded from below and the locality is 
guaranteed. In particular, for fx = 1 Ref. j^H found the condition 

O 2 

^->l-30e. (3.19) 
<r 

Even if this condition is not obeyed, numerical results can still reveal locality. 
In particular in the simulations at (3 = 6, 6.2, 6.4 the corresponding overlap 
operator can be made local [^j]]. The optimal locality is achieved by varying 
the parameter fx. The optimized values of fi at (3 — 6, 6.2, 6.4 were presented 
in Ref. [^j] and in Ref. jojj] the value of that guarantees a good locality at 
(3 = 5.85 was pointed out. 

The mass fx separates doublers from the "physical" eigenvalues of the 
Wilson operator, which are projected onto the left arc of the GW fermion 
spectrum. This imposes a restriction on fx to be between < fx < 2 in the 
free case. In the interacting case the interval shrinks from both sides. The 
separation can be done unambiguously only if the spectrum of the Wilson 
operator has a window with a negligible eigenvalue density. The latter cannot 
be achieved for the Wilson operator at strong (3. Therefore one rather has 
to work in a range of moderate values of (3. For illustration of this behavior 
we show in Figure 13.21 |9_4( two plots of the distribution of eigenvalues of the 
Wilson operator for quenched QCD configurations at f3 — 5.6 and (3 = 5.4 on 
the lattice of volume 4 4 . The left plot shows a spectrum with such a window 
while the right plot does not have it. Note that these distributions are done 
for quite small lattices and due to the boundary conditions the configurations 
are smother than they would be in a larger volume at the same (3. Hence the 
window in the spectrum opens up at somewhat lower values of (3 than it 
would be expected for large lattice volumes. 

Once the window was identified one can use the freedom of varying the 
value of fx within the window to optimize the locality properties further. 
This can be done by comparing the exponential decays of the couplings of 
the overlap operator at different fx. It turned out that the locality breaks 
down once (3 is below the value about 5.7 [95]. 

According to Ref. fTM it is optimal for locality to put fx — 1.4 at (3 — 6 
and from Ref. fi = 1.6 at (3 = 5.85. In Figure l3~3l we show the locality 
property of the overlap Wilson operator at (3 — 5.85, fx = 1.6. We take the 
"taxi driver metrics" \\x — y||taxi to measure the distance on the lattice 

4 

Ik -ylltoi = ^ \xn - Vn\ ■ (3.20) 
The locality is measured in terms of f(r) 
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Figure 3.2: Eigenvalues of the Wilson operator at (3 = 5.6 and (3 = 5.4 on the 
lattice 4 4 . We see that at (3 = 5.4 the spectrum does not develop a window. 



i c 

" + 
- + 

0.01 ~ 
0.0001 ~ 

<f(r)> [ 

1e-06 ~ 
1e-08 
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Figure 3.3: Locality of the overlap Wilson operator at /3 — 5.85, fi = 1.6 on 
the lattice of volume 12 3 x 24. The exponential decay is clearly pronounced. 
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f(r) = max 1 1 (D ov ) xy ip y \\ , where r = \ \x - x \ | tax i and ^ = <W, • 

(3.21) 

We see in Figure 13.31 that the decay is safely exponential and therefore the 
operator is local. 
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Numerical 
hypercube 



treatment of the 
fermions 



In this Chapter we discuss the simulation of the hypercube Dirac operator. 
First we consider the hypercube operator itself. We tune it to decrease the 
remaining additive mass renormalization and minimize the violation of the 
GWR, i.e. to improve its chiral behavior. Secondly we use the hypercube 
operator as a kernel for the overlap operator. Since the hypercube operator is 
a truncated fixed point operator we expect that the favorable properties such 
as locality become much better than for the overlap Wilson Dirac operator. 

Note that the hypercube operator would be perfect if we considered it for 
the free case on the lattice of size 3 4 . This provides us with an ansatz for the 
parameters (Ao, • • ■ , A4, k±, . . . , k 4 ) in Eq. ()2.143j) . We give the correspond- 
ing values in Table 14.11 On a larger lattice and for the interactive case we 
adopted the "minimal gauging" described in Subsection 2.5.4. As an ansatz 
for the couplings we take the values from the free case to start with. However 
this creates a large additive mass renormalization that was pointed out in 
Ref. |9fij |. In order to account for it we rescaled the parameters Aj and Kj as 
it was prescribed in Ref. |9fi| . 



where 1/u < v < 1. 

If we represent the hypercube operator as D HF [x, y, U] = 7 At p A1 (x, y, U) + 
X(x, y, U) then the rescaling of the couplings would imply a rescaling of the 



Ao 
Ai 
A 2 
A 3 
A 4 




Ao , 

M 3 A; 

u 4 A 



(4.1) 
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link variable t/ (U ( follows 

U XtlM -> uvU x>fl in p^{x,y,U) , 

U Xill -> uC/^ in X(x,y,U) . (4.2) 

The factor k is used to compensate the (mean) link suppression due to the 
gauge field and therefore reduce the mass renormalization. The factor v con- 
trols mostly the imaginary part of the eigenvalues and serves to move them 
closer to the GW circle which we described below Eq. (12.1301) . 

A further chiral improvement which was outlined in Ref. |S6J| is the use of 
the fat link. For a given configuration we substitute each link variable U Xtfl 
by the following expression 

(4.3) 

where the staples for the link U x ,p, are depicted on the Figure 14.11 The 
staples are calculated as a product of the three links and are given by: 
U XjV U x+Ptfl Ul + ^ v in positive direction and U x _ ^ v U x -p tfM U x -p+fi,u in negative 
direction. In the four dimensional space there are six staples for each link 
variable, hence the factor 6 in the definition. 

We applied a step by step procedure to tune couplings of the hypercube 
operator for quenched configurations generated at (3 = 5.85. We started from 
a 8 4 lattice and worked out a first approximation of a, u and v by plotting 
the eigenvalue spectrum of the hypercube operator Dhf and comparing it 
with the GW circle. In this way we obtained the following approximations 

a = 0.3 , u = 1.32 , v = 0.76 . (4.4) 

In Figure 14.21 we plot the lowest 80 eigenvalues of the hypercube operator 
for a typical configuration on 8 4 lattice at /3 — 5.85 for the mentioned pa- 
rameters. To improve the chiral properties of the hypercube operator further 
we decided to look at the condition number of the Q 2 operator resp. the 
polynomial degree of the overlap hypercube Dirac operator after projecting 
out 44 eigenvectors as a function of a, u and v for a given lattice volume 
L 3 x T. The accuracy of the Q 2 eigenvectors was taken 10~ 12 and the accu- 
racy for the polynomial approximation was of 10~ 15 . These accuracies were 
sufficient to obtain a good convergence of the Chebyshev polynomial ap- 
proximation. The parameter /i of an auxiliary mass in the definition of the 
operator Q = c7s(Dhf — A 4 ) c.f. Eq. (j3.3j) becomes redundant in this context 
since we can account for it by the proper combination of u and v. Therefore 
we set fj, = 1. We calculated the average value of the Q 2 condition number 
over 30 configurations — which we found to be sufficient for this purpose 



U x<fl ->• (1 - a)U Xtli + - ^ staples 
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- and built the gradient of it in the three dimensional parameter space of 
a, u and v. By going in the direction of the gradient, i.e. the direction of 
the steepest descent, we found the parameters optimized with respect to the 
condition number of Q 2 . This procedure was iterated until it converged. The 
optimal parameters are given in Table l4~2l 

Alternatively, one can also think of an optimization with respect to the 
locality properties. On the lattice 10 4 we found that if we increase the param- 
eter a then the locality properties of the corresponding overlap hypercube 
Dirac operator become worse. Therefore we decided to relax the optimal pa- 
rameter a for the condition number in favor of the locality of the resulting 
operator. We want to point out that the parameter a was decreased only 
mildly to provide the better locality property and still keeping the condition 
number at a good value. The suitable parameter a in this respect is a = 0.52 
on 10 4 lattice at (3 = 5.85. When we increased the lattice volume to 12 3 x 24 
we found out that the trend established on a smaller lattice is reproduced to 
a good extent. Therefore we can safely adopt the parameters from the study 
on the 10 4 lattice. This we first confirmed by the recalculating the gradient 
for the optimal parameters on the 12 3 x 24 lattice. Secondly it is supported 
by the behavior of the locality properties on 10 4 and 12 3 x 24 lattices. To 
illustrate this we show in Figure 14.31 the locality of the overlap hypercube 
Dirac operator at two different values of a together with the locality of the 
overlap Wilson Dirac operator. We see that the larger value of a, which cor- 
responds to a lower value of the condition number, shows a weaker slope and 
thus results in a worse locality. We therefore decide to stick in our simulations 
of the overlap hypercube Dirac operator to a = 0.52. We also observe that 
the locality behavior on the 10 4 lattice is reproduced almost identically on a 
larger 12 3 x 24 lattice. 

In Table FQ1 we give the condition numbers for the hypercube Q 2 operator 
for two optimized sets of parameters: with respect to the best condition 
number and the locality properties. We include the polynomial degree of 
the overlap hypercube operator after 44 eigenvectors being projected out at 
accuracy e over i ap = 10~ 15 . In Table l4~4l we add for comparison the condition 
number for the Wilson Q 2 operator. Note that the degree of the Chebyshev 
polynomial approximation is proportional to the square root of the condition 
number. 
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A 


1.8527205471652000 


Ai 


-0.0327938644185406 


A 2 


-0.0162118524304803 


A 3 


-0.0086184721387295 


A 4 


-0.0045483448675513 


«i 


0.0738626203107888 


^2 


0.0173136117866899 


^3 


0.0059685910389885 


K4 


0.0025632525141047 



Table 4.1: The parameters of the HF. 



^x-\-v,t& 



x-\-fi,i/ 



x ' x — // 

Ux,lj, 



Figure 4.1: A staple in the positive direction for the link variable U x 



+-Hr + + ^_ 



Figure 4.2: The lowest 80 eigenvalues of the hypercube operator Dhf with 
a = 0.3, u = 1.32, v = 0.76 on a 8 4 lattice, at ft = 5.85. For comparison we 
also show the corresponding GW circle. 
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V 


u 


V 


a 


8 4 


1.22 


0.76 


0.4 


10 4 


1.28 


0.96 


0.72 



Table 4.2: The parameters of the HF optimized with respect to the condition 
number at (5 — 5.85. 
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Figure 4.3: The locality of the overlap HF and the overlap Wilson fermions 
vs. the taxi driver distance r. The upper Figure shows the results on a 10 4 
lattice, and the lower Figure shows two results on 12 3 x 24 lattice along with 
10 4 lattice. We see the same exponential decay as we increase the lattice 
volume. 
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V 


u 


V 


a 


polynomial 


condition number 










degree 


m = 20 


m = 30 


m = 40 


8 4 


1.22 


0.76 


0.4 




25.2 ±0.6 


22.2 ±0.4 


20.2 ±0.3 


10 4 


1.28 


0.96 


0.52 


40 ± 1 


26.8 ±0.9 


22 ±0.5 


19.4 ±0.4 


10 4 


1.28 


0.96 


0.72 


39 ±1 


26 ± 1 


20.9 ±0.6 


18.2 ±0.3 


12 3 x 24 


1.28 


0.96 


0.52 


53 ± 1 


81 ±4 


52 ±3 


41 ± 1 



Table 4.3: The condition number for the hypercube Q 2 operator at fi = 1 
after m = 20, 30, 40 eigenvectors being projected out and the polynomial 
degree for the overlap hypercube Dirac operator after 44 eigenvectors being 
projected out. 



V 


polynomial 


condition number 




degree 


m = 20 


m = 30 


m = 40 


8 4 




417 ±10 


363 ±8 


330 ±4 


10 4 


172 ± 1 


636 ± 14 


529 ± 10 


467 ±7 


12 3 x 24 


247 ±4 


2100 ±200 


1260 ± 60 


996 ± 36 



Table 4.4: The condition number for the Wilson Q 2 operator at ji — 1.6 after 
m = 20, 30, 40 eigenvectors being projected out and the polynomial degree 
of the overlap Wilson Dirac operator after 44 eigenvectors being projected 
out. 
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Results on the pion dispersion 
relation 

The p-expansion can be used once the pion Compton wavelength is well 
accommodated inside the physical volume. This implies that in order to reach 
the chiral limit, i.e. very small quark masses, one has to simulate in extremely 
large lattice volumes. On the other hand the e-regime is defined when the 
pion Compton wavelength extends over the sizes of the physical volume. One 
uses the fact that the finite size effects appearing in this regime are controlled 
by the low energy constants of the infinite volume. In this Chapter we explore 
the product Lm n in the definition of the p- and e-regime. 

The pion mass in the infinite volume is extracted from a time exponential 
decay of the pseudo-scalar correlator a 3 t)P b (0, 0)). Here summa- 

tion over b is not taken and the pseudo-scalar density V b (x) was defined in 
Eq. (|2.62p . In the p-regime on a lattice of size L 3 x T this is modified to 

a 3 5> 6 (f,t)P 6 (0,0)> = l(0|^V° n)|2 e-^fcosh 

x n 

As the bare quark mass m q is decreased, the contribution of the zero modes to 
the (VV) rises as \v\ / (m q V) 2 , which is seen from the spectral representation 
of the pseudo-scalar correlator in Eqs. 1)2.65112.67)) . On the other hand \u\ ~ 
yV (we anticipate the formula (|7.5jl for the topological susceptibility). In 
total the contribution of the zero modes is reduced as the physical volume V is 
increased. In full QCD simulations the configurations with topological charge 
v would be suppressed compared to the quenched case due to the fermionic 
determinant by a factor m <T' • Therefore in our quenched simulations 
the sampling of the pseudo-scalar correlator has a systematic error hidden 
in the contribution of the zero modes. This error is increased for smaller m q 
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and decreased for larger V. As it was pointed out in Ref. |22|, the problem 
can be avoided by subtracting the scalar correlator (S b (x)S b (y)), 

S b (x) = $(x)iT b ip(x) , (5.2) 

which leads to the cancellation of the zero mode contributions. This hap- 
pens due to their exact chirality. The proposed method introduces also an 
additional scalar mass to the spectrum of the (VV — SS) correlator. It is, 
however, expected to be much larger than the pion mass and therefore it 
should be visible only among the excited states which are usually suppressed 
at least by a factor of exp (m n — -^excited)* - 

We used the Neuberger operator as a regularization for fermions to extract 
the pion mass in quenched QCD simulations on 12 3 x 24 lattice at (3 = 5.85. 
To calculate the pseudo-scalar correlation function the Dirac operator has to 
be inverted as in Eq. ()2.63|) at a bare quark mass m q which is introduced as 



£>ovK) = ( 1 - ^ ) Aw + m q . (5.3) 

We approximated the square root appearing in D ov with the Chebyshev poly- 
nomials to accuracy 1CT 15 . Descriptions of h igh l y op timized algorithms to 
invert the overlap operat or is given in Refs. |98l - ll0l| . In particular we use 
the multiple mass solver 102]. The quark masses were chosen to lie in the 



range where the p-regime is applicable so that the exponential decay of the 
pseudo-scalar correlator can be used. The lowest considered quark mass cor- 
responds to the cross-over of the p-regime to the e-regime. In Table 15.11 and 
in Figure l5~Tl the quark masses and the pion mas ses extracted from the VV 
as well as VV — SS correlators are reported [l03]. We see that the deviation 
between the two methods becomes larger as the quark mass is decreased. We 
also give the product m w L to quantify the the p-regime. In Figure IS~T1 we also 
show the chiral extrapolation by a linear fit of the square of the pion mass 
vs. the quark mass down to m q = 0. We see that the linear behavior is repro- 
duced well for the pion mass calculated with the (VV — SS) correlator. The 
intercept amounts to —0.002(6). For the lowest quark mass m q = 0.01, which 
is in the domain of the cross over to the e-regime, we found m n L = 1.68. This 
is quite small and therefore the finite size effects on the pion mass cannot be 
neglected. 



To fix the physical scale we used the Sommer scale r n, see Refs. |104l .ll 
It was computed for the considered range of (3 in Ref. |l06l | 
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Figure 5.1: The square of the pion mass nv\ vs. the bare quark mass computed 
both for (W — SS) and (VV) correlators. The data for (VV — SS) can be 
linearly extrapolated to m q a = with an intercept compatible with m^a 2 = 
0. 



65 



Chapter 5 Results on the pion dispersion relation 



m q a 


m' w a 


v-s 
m n a 




0.01 


0.212(9) 


0.140(20) 


1.68 


0.02 


0.237(7) 


0.196(14) 


2.35 


0.04 


0.299(5) 


0.280(10) 


3.36 


0.06 


0.355(4) 


0.346(8) 


4.15 


0.08 


0.405(4) 


0.401(7) 


4.81 


0.10 


0.450(4) 


0.451(6) 


5.41 



Table 5.1: The pion mass calculated from 140 configurations for the 
{VP) and for the (VV — SS) correlators. The products m n L are also shown. 
At m q a = 0.01 we obtain the pion mass m^~ s rs 225 MeV. 
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Eigenvalue distributions of the 
overlap Dirac operator 



Since the individual distributions of the eigenvalues of the Dirac operator 
were only conjectured in the framework of the xRMT, it is of interest to 
compare these predictions to the corresponding eigenvalue distributions ob- 
tained from QCD lattice simulations with chiral fermions. 

The Wilson Dirac operator seems hardly suitable in this respect; because 
of the mass renormalization it would be highly problematic to identify the 
zero modes and the relative values of the remaining modes. In addition, 
quenched simulations at small quark masses are plagued by the occurrence 
of accidentally extremely small unphysical eigenvalues which is known as 
exceptional configurations. 

Staggered fermions do not suffer from additive mass renormalization be- 
cause they have an exact remnant chiral symmetry £7(1) <8> U(l) . Indeed 
such comparisons exist. The staggered fermion spectrum agrees well with 
the xRMT prediction in the sector v — 0. However, it turns out that all 
other sectors, v ^ 0, yield t he same distributions, in particular the same 
histogram for Pi\z) in QCD (l07t Il08| ] . It seems that staggered fermions are 
generally insensitive to the topology, and therefore not adequate in their con- 
ventional formulation for this purpose, at least for moderate and strong gauge 
coupling. However, latest studies with improved staggered fermions, where 
the mixing of the pseudo-fl avors is suppressed, showed that the sensitivity 
can be recovered |l03 - llll| . There is also an indication that the sensitivity 
of the conventional staggered fermions becom es b etter as we approach the 
continuum limit. This was suggested in Ref. 112l | where a study with the 
Schwinger model was presented. 

The lattice studies of the distrib ution s of individual eigenvalues with the 
GW fermions were pioneered in Ref. 113] where the eigenvalues of the overlap 
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operator and the fixed point operator were considered for the Schwinger 
model. In both cases the results agreed with the xRMT formulae. 

Su ch s tudies were also performed with the overlap operator in 4d 



QED |ll4j and again the predictions were successful within the statistical 
errors. 

F inall y QCD was considered, but initially just on tiny lattices of size 4 . 
Ref. 1151 ] used the overlap operator at strong coupling of f3 = 5.1 which in our 



opinion corresponds to too strong coupling regime where the ov erlap formula 
cannot be safely applied, as we pointed out in Chapter 3. Ref. O applied 



a truncated fixed point action and obtained a decent agreement in a volume 
of V = (1.2 fm) 4 . However, when the physical lattice spacing is decreased 
so that the volume shrinks to (0.88 fm) 4 , the leading non-zero eigenvalue 
distributions of the different topological sectors are on top of each other, in 
contrast to the xRMT prediction. 



6.1 Microscopic regime 

In this Section we discuss our results for the individual eigenvalue distribu- 
tions ofthe overlap Wilson Dirac operator. Most of them were published in 



Ref. 117 



Due to the computational time required by the overlap operator we had to 
use the quenched approximation, as it was also the case in all previous studies 
mentioned above. We generated quenched SU(3) gauge fields with the Wilson 
gauge action ()2.110|) with periodic boundary conditions. We used moderate 
gauge couplings, (3 > 5.85 and chose the mass fi in the overlap formula in 
a way to render the operator optimally local. In particular we took u. — 1.4 
for (3 = 6 and fi = 1.6 for (3 = 5.85. We approximated the inverse square 
root by Chebyshev polynomials to an accuracy of e over i ap = 10~ 1 2 for latti ces 
of volume 8 4 , 10 4 , 12 4 . Then we used the PARPACK routines [ill El to 



evaluate up to 100 eigenvalues of the overlap Wilson Dirac operator with 
the least absolute value. The index for each configuration was identified by 
computing eigenvalues of the operator 75 for the leading modes. They came 
out to be ±1 to a very high accuracy for the zero modes. The corresponding 
zero eigenvalues on the 10 4 lattice at (3 = 5.85 were of order 10~ 10 to 10~ 6 , and 
the first non-zero eigenvalues were in the range about 10~ 3 to 10~ 2 . So there 
was a large gap between non-zero and zero eigenvalues, which allowed us to 
safely identify the zero modes. Then the index was read off from the chirality 
of the corresponding eigenvectors. To illustrate this we show in Table 16.11 an 
example of the eigenvalues and the chirality of the corresponding eigenvectors 
of the Neuberger operator in topological sectors v = and —2. Of course 
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EV of DlD ov 


EV of 7s 


example for charge v — 


6.23013e - 3 


0.3608201 


6.23014e - 3 


-0.3608198 


9.90174e - 3 


-0.7622260 


9.90177e - 3 


0.7622260 


2.69086e - 2 


-0.7528350 


example for charge v = — 2 


6.03761e- 16 


-1.0000000 


9.43764e- 11 


-0.9999999 


3.68203e - 3 


-0.9240370 


3.73623e - 3 


0.9205659 


7.04176e - 3 


-0.7016518 



Table 6.1: Typical examples for the lowest eigenvalues of Dl v D ov on a 10 3 x24 
lattice at (3 = 6. 



lattice 
size 


P 




physical 
volume 


numbe 
v = 


r of confi 
\u\ = 1 


gurations 
\u\ = 2 


8 4 


5.85 


1.6 


(0.98 fm) 4 


80 


63 


28 


12 4 


6 


1.4 


(1.12 fm) 4 


44 


70 


24 


10 4 


5.85 


1.6 


(1.23 fm) 4 


74 


112 


84 


16 3 x 32 


6 


1.4 


(1.49 fm) 3 x 2.98 fm 


49 


87 


63 



Table 6.2: The statistics of our simulations on four lattice sizes. 



the non-zero eigenvalues occur in complex conjugate pairs, hence we will 
consider only one sign for the imaginary part below. Our statistics in various 
topological sectors is collected in Table 16.21 To fix the physical scale we again 
used the Sommer scale introduced in Chapter 5. 

6.1.1 Distributions of individual eigenvalues 

In this Subsection we discuss the probability distributions of the individual 
lowest eigenvalues of the overlap Wilson Dirac operator. 

In order to relate the eigenvalues found on the Ginsparg- Wilson circle to 
the continuum eigenvalues A n , we map the circle stereographically onto the 
imaginary axis. Requiring f(z) = z + 0(z 2 ) and /(2/i) = oo singles out the 
Mobius transform 
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This mapping has been suggested before in Refs. |120L Il21| . for instance, in 
connection with the Leutwyler-Smilga sum rules [120J. 

The results for the distribution of the lowest eigenvalue of the overlap 
Wilson Dirac operator can be represented by a histogram, which then could 
be compared to Figure l2~Tl However, such a picture depends on the arbitrary 
choice of the bin size in the histogram. This can be avoided by calculating the 
cumulative density, which sums up all the entries up to the considered value 
of z — it is normalized so that the full number of entries corresponds to the 
cumulative density 1. The resulting cumulative densities can be compared to 
the analytical curves of the cumulative density of the xRMT prediction for 
the kth eigenvalue in the topological sector u, 



u<>, , _ Jo P { k \ z ') dz> 

S™pt\z')dz> 



The chiral condensate enters these analytical formulae as the only free param- 
eter. Therefore, applying a one parameter fit to the data in various topological 
sectors it is possible to obtain an estimate for its bare quenched value. 

We note in passing that the absolute values of the eigenvalues of the 
operator 7s-D ov are equal to absolute values of the a/ Dt v D ov . Moreover the 
stereographically projected eigenvalues of the ope rator D m can be related to 
the Hermitian eigenvalues by the simple formula 



10(1 |l22 



A, 



1 - V(2 M ) 



|Af| 



where (6.3) 



1 - (—\ H ) 2 

D ov ipj = iXjipj , 75-D ov 0j = i\f4>j ■ (6.4) 

We start the discussion of the cumulative distributions for our small- 
est physical volume of (0.98 fm) 4 in Figure I6~T1 Here the lines represent the 
predictions by xRMT for the lowest eigenvalue of the Dirac operator in the 
topological sectors \v\ — 0, 1 and 2, and the symbols the cumulative distribu- 
tion of the overlap Wilson Dirac operator. One can see a clear disagreement 
with the xRMT curves. In particular the sensitivity of the data to the topo- 
logical charge is lost to a large extent. Notice that for very small values of 
z the data are closer to the analytical predictions. However, this does not 
allow for a reasonable fitting procedure to extract the chiral condensate S. 

Now we consider a larger physical volume of (1.12 fm) 4 . It corresponds 
in our simulations to a 12 4 lattice volume at (3 — 6. The plot is depicted in 
Figure l6~2l One can see that the agreement starts to set in. We interpret this 
picture as a satisfactory agreement with the predictions by the xRMT. In 
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(3 = 5.85 
L = 8 

M = 2 



2 4 6 8 10 

z 

Figure 6.1: The cumulative distribution of the lowest non-zero eigenvalue of 
the overlap Wilson Dirac operator on 8 4 lattice, at (3 = 5.85 in the topological 
sectors \v\ — 0, 1, 2 vs. xRMT predictions (lines). No agreement is found. 




2 4 6 8 10 

z 



Figure 6.2: The cumulative distribution of the lowest non-zero eigenvalue of 
the overlap Wilson Dirac operator on 12 4 lattice, at (3 = 6 in the topological 
sectors \u\ = 0, 1, 2 vs. xRMT predictions at S = (256 Me V) 3 (lines). A 
reasonable agreement is observed. 
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o 



L= 10 



P = 5.85 



M = 2 



o 



2 



4 



6 



8 



10 



Z 



Figure 6.3: The cumulative distribution of the lowest non-zero eigenvalue 
of the overlap Wilson Dirac operator on 10 4 lattice, at /3 — 5.85 in the 
topological sectors \u\ = 0, 1, 2 vs. %RMT predictions at £ = (253 MeV) 3 
(lines). A good agreement is observed. 

particular, we do see the effect that the peak of the density, resp. the interval 
of steepest ascent of the cumulative density, moves to larger values of z for 
increasing topological charge \u\. The plot is shown for the optimal value of 
the chiral condensate £ = (256 MeV) 3 . 

Next we studied a 10 4 lattice at /3 = 5.85, which corresponds to a some- 
what larger physical volume of V = (1.23 fm) 4 . Again the results for the 
leading non-zero eigenvalue are shown for the topological sectors \v\ = 0, 
1 and 2 in Figure 16.31 We interpret the results to be in a good agreement 
with the xRMT predictions. The optimal value of the chiral condensate is 
modified to £ = (253MeV) 3 . 

Finally we looked at the lattice size 16 3 x 32, at /3 — 6 which corresponds 
to a yet larger volume of (1.49 fm) 3 x 2.98 fm. In Figure lfi~H we show the 
cumulative densities of the lowest eigenvalue in topological sectors \u\ = 0, 1 
and 2 for £ = (286 MeV) 3 . We interpret it again as a good agreement. 

yRM T predicts a (logarithmically) increasing E(V) due to quench- 
ing 123]. A conclusive verification of this behavior (beyond possible lattice 
artifacts) would require further simulations. 

In Figure 16.51 we show the cumulative distributions of the next-to-lowest 
eigenvalue in the topological sectors \v\ — 0, 1 and 2. The data we take 
from simulations on 10 4 and 12 4 lattices. Here we do not find a convincing 
agreement with the %RMT predictions; note that the relevant values of z are 
larger compared to Figure IB~2l Also the second non-zero eigenvalue moves to 
larger values of z if \u\ increases. 

To quantify the quality of the fit we used the Kolmogorov-Smirnov 
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2 4 6 8 10 



Figure 6.4: The cumulative distribution of the lowest non-zero eigenvalue 
of the overlap Wilson Dirac operator on 16 3 x 32 lattice, at (3 = 6 in the 
topological sectors \u\ = 0, 1, 2 vs. xRMT predictions at S = (286MeV) 3 
(lines). A good agreement is observed. 




Figure 6.5: The results for the cumulative distribution of the next to lowest 
eigenvalue in topological sectors \u\ = 0, 1, 2 for two lattice volumes 10 4 and 
12 4 vs. the xRMT predictions (lines). No convincing agreement is observed. 
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lattice 
size 




optimal 
S 


CO 

v = 


tifidence 1 
\v\ = 1 


evel 
\v\ = 2 


12 4 


6 


(256MeV) 3 


0.003 


0.73 


0.79 


10 4 


5.85 


(253MeV) 3 


0.03 


0.48 


0.10 


16 3 x 32 


6 


(286MeV) 3 


0.7 


0.8 


0.12 



Table 6.3: The results for the chiral condensate and the statistical confidence 
level according to the Kolmogorov-Smirnov test. 



test [92j. This test provides a confidence level that the given data set is 
drawn from the probability distribution of interest. In Table l6~3*l we collect 
the confidence levels and the optimal values for £ on the lattice volumes 
where we found a reasonable agreement with the xRMT. We see that the 
box length of 1.2 fm constitutes the lower bound of applicability of the xRMT 
predictions. We do not observe a full agreement. However, one can argue that 
the separation of the topologies for the lowest eigenvalue has set in to a large 
extent, and the box length is sufficient to obtain a good estimate for the 
chiral condensate. 

The observation that there is a minimal physical box length where the 
xRMT predictions set in can be interpreted in terms of the so-called Tho uless 



energy, a term that originates from mesoscopic physics (see the Review 124 
on more details). The xRMT predictions become valid once we ensure that 
the simulations are performed in the e-regime of the x?T. This can be esti- 
mated if we require the potential term in the Lagrangian of the %PT f)2.37j) to 
dominate over the kinetic term. Applying this argument to the characteristic 
scale of Dirac eigenvalues we arrive at 

F l 

A <C Axhouless — ) (6.5) 

or converting this to the dimensionless variable: ^Thouless = F^y/V. The above 
relation can also be obtained from the lower bound in the relation 1)2.3.2)1 
by taking its square and recalling the formula for the pion mass through 
the quark mass, see Subsection 2.3.1. The determination of the value of z 
where the xRMT predictions fail can be considered as an empirical value for 
the Thouless energy. Thus the Thouless energy sets a limit on the value of 
eigenvalues where the xRMT applies. On the other hand one has to ensure 
that the physical volume is large enough to be in the confinement phase 
where all assumptions for the %PT are made. In addition one has to keep 
in mind the relation ()2.43)) ensuring that the %PT expansion is defined. The 
latter imposes on the physical extent the requirement L ^> l/(2F 7r ) ~ 1 fm. 



74 



Chapter 6 Eigenvalue distributions of the overlap Dirac operator 



To check the condition ()6.5j) one can perform the same eigenvalue study at 
rather large fixed physical volume and different lattice spacings a. Due to the 
renormalization of and S this will lead to different Thouless energies and 
therefore it should be seen in the behavior of the cor respo nding cumulative 
distributions. In fact such a study was done in Ref. 
confirmation of this behavior was found. 



125l |. but no ultimate 



6.1.2 Spectral density 

In this Subsection we discuss the microscopic spectral density of the over- 
lap Wilson Dirac operator. We still separate topological sectors, but now we 
consider the distributions of all computed eigenvalues. We show the eigen- 
value distributions for topological sector \v\ — 1 and two lattice volumes, 10 4 
and 12 4 , in Figure xRMT predicts an oscillating behavior ()2.98j) 
which is also plotted for comparison. However, here we omit the zero eigen- 
values, hence the term with the delta function drops out. We find a good 
agreement roughly up to the second peak, as expected from Figures 16. 2| 16.31 
and 16.51 Then we are leaving the microscopic regime and turn to the bulk. 
The Thouless energy appears in this context as a value of z where we leave 
the microscopic regime. 

If the oscillation is averaged to a plateau, its height agrees well with the 
eigenvalue density at zero according to the Banks-Casher relation, p(z = 
0) =p(A = 0)/EV = 1/tt. 



6.2 Bulk eigenvalues 

In this Section we discuss the eigenvalues in the bulk of the Dirac spectrum. 
It is expected that the spectral density in the deep bulk domain has a A 3 
behavior. Since we were able to compute as many as 50 conjugate pairs of 
eigenvalues of the overlap Wilson Dirac operator we can analyze the spectral 
density in the bulk. The largest eigenvalue corresponds to an energy of about 
800 MeV. However, it is difficult to argue whether the considered domain 
is sufficiently large to accommodate the A 3 behavior. In Figure 16.71 at the 
bottom we plot the eigenvalue spectral density for the topological sector 
v = 1 on the lattice 12 4 against the A 3 curve, which we take to be a + bz 3 . 
The data follow amazingly well the predicted behavior. 

Now we want to make a remark on the Banks-Casher relation. If we were 
not able to obtain information from the microscopic regime and we had only 
the eigenvalues in the bulk, then the chiral condensate which would follow 
from the approximation of the eigenvalue density would be overestimated. 
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Figure 6.6: The results on the microscopic spectral density p s (z) in topo- 
logical sectors \u\ = 1 for two lattice volumes 10 4 and 12 4 vs. the xRMT 
predictions (lines). 
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Figure 6.7: The upper plot shows the spectral density p s (z) in the bulk vs. 
the microscopic regime. One line corresponds to a + bz 3 and another to the 
microscopic spectral density. The plot at the bottom shows the spectral den- 
sity in the bulk. The cross-over from the microscopic regime to the bulk is 
nicely seen. 



This can be seen from the plateau extrapolated to low values of z on the 
Figure 16.71 on the top. For illustration we plot two curves for the bulk and 
the microscopic regime. 
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6.2.1 The unfolded spectrum 

We now turn to another way of comparing our lattice data to a conjecture 
from %RMT. This additional evaluation allows us to take all our non-zero 
eigenvalues into account (for one sign of the imaginary part, i.e. up to 50 for 
each configuration). We build from all the eigenvalues of all configur ation s 
the hierarchically unfolded distribution as described for instance in Ref . |ll5| . 
For a co mpr ehensive review and other ways of unfolding we refer for instance 
to Ref. [l24 |. 

To construct the hierarchically unfolded spectrum we first numerate all 
available non-zero eigenvalues with positive imaginary part in each configu- 
ration in ascending order, given by the angle in the Ginsparg- Wilson circle. 
Then we put all these eigenvalues from all configurations together and nu- 
merate them again in ascending order. Now we consider pairs of eigenvalues 
from the same configuration, which follow immediately one after the other 
in the original numeration. If they differ in the global numeration by k, 
then k/N con { is the hierarchically unfolded level spacing, where iV con f is the 
number of configurations involved. In Figure IH7H1 we show our results for the 
corresponding level spacing distribution in V = 12 4 , (3 = 6, and in V = 8 4 , 
(3 = 5.85. The histograms are compared to RMT conjectures for different 
groups and representations, and we can confirm a good agreement with the 
Wigner distribution predicted for SU(3) in the fundamental representation. 
We conclude that this property holds over a large range of volumes. Since we 
include all our eigenvalues here, our statistics is much larger than in the plots 
for the microscopic regime. Note, however, that this analysis is not sensitive 
to the topology any more. There is no physics information such as values 
for the condensate or the pion decay constant, because one does not keep 
track of a physical scale. The hierarchically unfolded distribution only tests 
the correct symmetry group that is to be taken in RMT. Nevertheless, the 
hierarchically unfolded distribution provides a non-trivial and successful test 
for xRMT. 
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Figure 6.8: The hierarchically unfolded level spacing distributions. For com- 
parison we show the curves for SU(2) and SU(3) fundamental representation 
and for SU (3) adjoint representation. 
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7.1 Motivation 

The topological susceptibility for a Yang-Mills gauge theory is defined by 

Xt = J d 4 xd 4 y (q{x)q{y)) Y M , (7.1) 

where q(x) is a density of the topological charge Q = J d xq(x) and the 
subscript YM means that the average is normalized by the partition func- 
tion. We refer to quenched QCD, although the topological susceptibility also 
exists in the dynamical case. The main reason to look at the quen ched topo- 



logical susceptibility comes from the Witten-Veneziano formula [126, 127 . 
It is motivated by the inability of the constituent quark model to describe 
a large mass of the rj' meson. The Witten-Veneziano formula predicts that 
at the leading order in the Nf/N c expansion, with Nf being the number of 
flavors and N c the number of colors, the contribution to the rj mass due to 
the Ua(1) anomaly is given by 

m% = ^fxt ■ (7.2) 

In the absence of the anomaly one could argue that Ua(1) symmetry 
would be broken spontaneously leading to a massless Goldstone boson for 
massless quarks. However, a particle with the corresponding quantum num- 
bers is not found, which is attributed to the existence of the Ua{X) anomaly. 

To interpret Eq. (j7.1j) one has to find a properly normalized topologi- 
cal density q(x) and to subtract from the quantity q(x)q(0) an appropriate 
contact term to define it properly in the continuum limit. 
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In the case of Ginsparg- Wilson fermions the topological density can be 
determined by 

q x = -Ltt[ 15 D[UU x ] , (7-3) 
2a 6 

where the trace is taken over color and spinor indices of a Ginsparg- Wilson 
operator D. It can be shown that this definition leads to a pr operly defined 
topolog ical charge and to the Atiyah-Singer index theorem (^J Il28l . E^29j . 
In Ref. [l30] it was also demonstrated that this discretization naturally enters 
the Witten-Veneziano formula if we consider Ginsparg- Wilson regularization 
for the fermions 



2N 



f 



F 2 



hm a 

V-*oo V 



77 y^(gygx)YM 



(7.4) 



x.y 



Note that due to the chiral symmetry of the fermions no additional sub- 
traction is necessary |l3lj . Using Eq. ()7.3|) we obtain for the topological 
susceptibility 

((N + - iV_) 2 ) YM 



Xt 



lim 

V- 
a 



-0 



V 



(7.5) 



where N + and AL is the number of zero modes of the Ginsparg- Wilson 
fermions with positive and negative chirality, respectively. This formula will 
be a starting point for the discussion of the results which we obtained using 
the overlap operator. 

The study of the continuum limit o f the topological s usceptibility with 

the GW fermions was pioneered in Refs. Ea, E33] it is argued that 

the continuum extra polat ed value of the topological susceptibility is about 

1331 ] . where r = 0.5 fm. 



Xt r\ = 0.059 ± 0.003 



7.2 Results with overlap fermions 

We simulated quenched QCD on I 3 x T lattices with periodic boundary con- 
ditions and used the overlap Wilson Dirac operator, as well as the overlap 
hypercube Dirac operator, to determine the index of the gauge configura- 
tions. The configurations were simulated using the usual Monte Carlo pro- 
cedure where we employed the heat bath and the overrerlaxation method 
as described in Subsection 3.1. We assume the configurations to be statisti- 
cally independent, which we checked by measuring the autocorrelation time 
of the topological charge. It turned out to be below 0.5, which characterizes 
a completely decorrelated sample. 
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Histogram of indices on a 16 x 32 lattice at p = 6 with 401 configurations 




Figure 7.1: The charge histogram and the corresponding Gaussian distribu- 
tion with the same (is 2 ) and norm. 



We considered several lattice volumes and in Table 17.11 we collect the 
topological susceptibility for the Neuberger operator as well as for the overlap 
hypercube operator. 

One expects a Gaussian charge distribution which is confirmed by our 
data, see Figure 17.11 In order to render the values of the topological suscep- 
tibility stable against the fluctuations due to a few entries of high charges 
(\v\ > 8) on the lattice volumes 12 3 x 24 and 16 3 x 32 we applied the following 
method: we counted the number of entries q v in the interval [—is — 0.5, is + 0.5] 
and normalized them to the total number of entries. Then we employed the 
formula 



Qv 



u 



/27T 



exp(-y 2 /2)dy 



0.5 



(7.6) 



to extract the value of (is 2 ) and thereby the topological susceptibility. We 
calculated it for \is\ = 1,2, ... , 5. Out of these values we computed a mean 
value and a standard deviation. In Table 17.21 the topological susceptibility 
calculated with this method is shown. 

There are two issues when considering the topological susceptibility calcu- 
lated on the lattice. First the finite size effects dependence of the topological 
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susceptibility, which is governed by exp (— LM g i ue b a ii), where M g i ue baii ~ 1-5 
GeV is the mass of the lightest glueball. Thus for lattices of L > 1 fm the 
finite size effects on the topological susceptibility are expected to be far be- 
low our statistical errors. We see from Table 17711 that at (3 — 5.85 and lattice 
volumes V = 10 4 and 12 3 x 24 the corresponding values of the topological 
susceptibility coincide within the error bars. Second, at finite lattice spac- 
ing a the topological susceptibility suffers from discretization effects starting 
from 0(a 2 ). 

In Figure 17.21 we plot the topological susceptibility calculated with the 
Neuberger operator and the overlap hypercube operator using Eq. ()7.6|) vs. 
the square of the lattice spacing. We fixed the physical volume to (1.48 fm) 3 x 
2.96 fm. We show the results for the Neuberger operator for 12 3 x 24 volume 
at (3 — 5.85 and 16 3 x 32 at (3 = 6. The result for the overlap hypercube 
fermions is given for the volume 12 3 x 24 at (3 = 5.8 5. Fo r comparison we 
also show the continuum extrapolated value from Ref . |l33| . We see that the 
value of the topological susceptibility calculated with the overlap hypercube 
Dirac operator is closer to the continuum value than the corresponding value 
of the Neuberger operator. This result hints at an improved scaling behavior 
of the overlap HF, though further statistics and observables are needed to 
certify this trend. 
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Topological susceptibility in V = (1 .48 fm) 3 x 2.96 fm 
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Figure 7.2: The topological susceptibility Xt r o> calculated with the Neuberger 
operator and the overlap hypercube Dirac operator, vs. the square of the 
lattice spacing. 



operator 


lattice 
size 


P 


n 


physical 
volume 




iVconf 


Neuberger 


12 4 


6 


1.4 


(1.12 fm) 4 


0.063 ±0.008 


143 


Neuberger 


10 4 


5.85 


1.6 


(1.23 fm) 4 


0.078 ±0.006 


308 


Neuberger 


12 3 x 24 


5.85 


1.6 


(1.48 fm) 3 x 2.96 fm 


0.068 ± 0.0048 


400 


Neuberger 


16 3 x 32 


6 


1.4 


(1.49 fm) 3 x 2.98 fm 


0.073 ± 0.0053 


401 


overlap HF 


12 3 x 24 


5.85 


1 


(1.48 fm) 3 x 2.96 fm 


0.066 ± 0.0043 


391 



Table 7.1: The topological susceptibility calculated with the overlap hyper- 
cube Dirac operator and with the Neuberger operator using the standard 
deviation. 



operator 


lattice 
volume 


P 


fi 


physical 
volume 




-^conf 


Neuberger 


12 3 x 24 


5.85 


1.6 


(1.48 fm) 3 x 2.96 fm 


0.069 ±0.0013 


400 


Neuberger 


16 3 x 32 


6 


1.4 


(1.49 fm) 3 x 2.98 fm 


0.072 ± 0.0054 


401 


overlap HF 


12 3 x 24 


5.85 


1 


(1.48 fm) 3 x 2.96 fm 


0.065 ±0.0027 


391 



Table 7.2: The topological susceptibility calculated with the overlap hyper- 
cube Dirac operator and the Neuberger operator using Eq. (|7.6|) 
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8.1 The axial- vector correlator 

In this Section we discuss results on the axial-vector correlator. Most of these 



results were published in Ref. |134j. We simulated quenched QCD with the 
Wilson gauge action on the lattices of size 10 3 x 24 and 12 4 at (3 = 6 with 
periodic boundary conditions. As a regularization scheme for the fermions 
we chose the overlap Wilson Dirac operator. We calculated the axial-vector 
correlator according to formula ([2.54)1 . To this end we needed to invert the 
overlap Wilson Dirac operator at a small quark mass m q . We approximated 
the square root appearing in D ov with Chebyshev polynomials to accuracy 
1CT 16 . The statistics of the configurations for various topolo gical secto rs is 
given in Table IH7T1 To cony Dute the index we used PARPACK [23, Ojjj and 



Ritz functional methods |135i Il3fi| . We kept the quark masses light enough 
in order to stay in the e-regime. To quantify this condition, z q = m q Y>V must 
be well below 1. In particular the considered quark masses were chosen in a 
way that z q < 0.34. We determined the parameter F n by making a fit with an 
additive constant in Eq. (|2.61j) at r = 0.5. Thus it can easily be computed. 
The chiral condensate E, on the other hand, is proportional to the curvature 
of the parabola. However, as the value of E varies in the range of interest, the 
corresponding parabola changes only very little which makes the extraction 
of E extremely difficult. 

We find that the lattice size 10 3 x 24 at (3 = 6 is too small to accommodate 
the predictions by the ;^PT. In Figure 18.11 we show the data for the axial- 
vector correlator at quark mass m q = 21.3 MeV, lattice volume V = 10 3 x 24, 
at (3 = 6 for topological sectors \v\ = 1,2 along with the corresponding ana- 
lytical predictions. We plot the analytical curves for two very distinct values 
of the chiral condensate, E = and E = (250 MeV) 3 . We see that the re- 
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Figure 8.1: The axial- vector correlator computed at quark mass 21.3 Me V 
for a lattice volume 10 3 x 24 at (3 = 6. The left plot shows the axial- vector 
correlator for topological charge v = and the right plot for \v\ = 1. The 
result at v = suffers from large error bars. In the sector \u\ = 1 we see a 
clear failure to reproduce the qxPT predictions (curves). 



suiting curves lie too close to each other to allow for E to be evaluated. The 
fact that those two curves are almost indistinguishable is valid for topological 
sectors \u\ ^ 0. It is due to the last term in the expression (|2.61j) . However in 
the topological sector the sensitivity is restored. At \u\ = 1 we see that the 
analytical prediction for this lattice volume does not model the data. The 
failure of the xPT to describe the result on volumes with box length less 
than 1.2 fm has been observed before in Section 6.1 where the distributions 
of eigenvalues of the overlap Wilson Dirac operator were confronted with the 
xRMT predictions. Thus this conclusion is rather general and sets a consis- 
tent constraint in the lattice simulations where the results can be confronted 
with xPT. 

Next we consider the fit by the analytical prediction of the %PT to the 
simulation data on the volume 12 4 and topological sector \u\ = 1. Our xRMT 
study of Dirac eigenvalues in Chapter 6 suggests that this volume should be 
at least at the threshold where the predictions by xPT become valid. The 
result is shown in Figure 18.21 (on the left) where we fit Eq. (j2.61|) to the 
simulation data, omitting the points located far from the center which are 
strongly affected by the contributions of the excited states. We can compute 
the additive constant that corresponds to F n quite accurately. Figure 18.21 
(on the right) shows the fitted value of F n against the extension of the time 
interval where the fit was done. We find a decent plateau which suggests a 
value F n = 87 ±4 MeV for the quenched, bare pion decay constant. Again the 
chiral condensate cannot be extracted from these data. This we demonstrate 
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Figure 8.2: The fit of the quenched QCD simulation data for the axial-vector 
correlator vs. the analytical predictions by xPT (left plot). The curves rep- 
resent the qx?T predictions for E = and E = (250 MeV) 3 . The extracted 
value of F n vs. number of points tf for the fit range around t/T = 0.5 (is 
shown on the right plot). The lattice volume is 12 4 , (3 — 6, \u\ = 1 and 
m q = 21.3 MeV. 




in Figure l8~2"l fon the left) where we plot the analytical prediction for two cases 
E = and E = (250 MeV) 3 . We see that the resulting curves reside within 
the error bars and hence make the determination of E impossible. Finally we 
look at the simulation data on 12 4 lattice at (3 = 6 for the topological sector 
\u\ =2. To check whether the data agree with the analytical predictions we 
varied F n in the range that we found from the fit in the topological sector 
\u\ = 1. The resulting shaded area is shown in Figure I5~3l We conclude that 
the data for the topological sectors \v\ — 1 and 2 agree within the errors. 
Note that the agreement at \v\ — 2 is, however, somewhat worse than for 
the topological sector \v\ = 1 (the data are more flat). This is, in fact, in 
agreement with the xRMT study of the eigenvalues of the Dirac operator, 
where we found that the agreement with the ^RMT is expected to be worse 
as we increase the topological sector \u\. 

We obtained of course the bare values for F n which are still subject to 
the renor malization . The latter is expected to increase the bare values sub- 
stantially. |loill37| . 
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Figure 8.3: The quenched QCD simulation data for the axial- vector correlator 
vs. the analytical predictions with the parameter F n in the range that was 
found by fitting in the topological sector \u\ = 1 (shaded area). The lattice 
volume is V = 12 3 , at /3 = 6, for topological charge \u\ = 2. The quark mass 
is m q = 21.3 MeV. 

8.2 Subtleties of numerical simulations in the 
e— regime 

In this Section we address the practical side of the simulations in the e- 
regime. 

In the e-regime the sector of topological charge v — may appeal as a 
good candidate where the chiral condensate can be extracted since the ex- 
pression (J2.61)) does not show such an insensitivity to £ as it is the case for 
non-zero topological sectors. Unfortunately, this sector is not favorable for 
the numerical simulations because of the frequently occurring small eigen- 
values. As we saw in Section 2.4 and in Chapter 6, where we studied the 
probability distribution of the lowest eigenvalues, it is the topological sector 
v = where we observe a high probability of very small eigenvalues Ai. As we 
increase the topological charge, the corresponding probability is suppressed, 
see for instance Figures and ED This behavior is reflected in the relatively 
large error bars for the axial-vector correlator in the topological sector 0, il- 
lustrated in Figure 18.11 We can demonstrate this explicitly by showing how 
the axial- vector correlator is sampled for the topological sectors and 1. In 
Figure E31 the Monte Carlo histories of A (t)A (0)\ u=0 and A (t)A (0)\ M=1 
are shown. In particular we see strong spikes in the history of the axial-vector 
correlator for the neutral topological sector. We checked the configurations 
where they show up. They correspond to cases of small eigenvalues of the 
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Figure 8.4: Monte Carlo histories of the axial- vector correlator for topological 
sector v = (left plot) and topological sector \u\ = 1 (right plot), at V = 
10 3 x 24, p = 6. 
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Figure 8.5: The height of the spikes vs. the bare quark mass. 



overlap Wilson Dirac operator. The heights of the spikes become larger for 
smaller quark masses. This is shown in Figure 18.51 The height of the spikes 
rises monotonously in r^-r , where Ai is the lowest eigenvalue of the Dirac oper- 
ator. This also holds for the spikes in the non-trivial sectors, where, however, 
the occurrence of such very small modes is suppressed. In practice one has 
to trade off the quark mass in a way first to ensure that it is small enough to 
be in the e-regime and secondly to have it large enough to damp down the 
values of the spikes. 

Now we want to trace the effect of the spikes on the statistics necessary 
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Figure 8.6: The standard deviation of £^/ n for topological sectors \v\ = 0, 1 
and 2 vs. pseudo statistics. The stabilization is observed much earlier for the 
topological sectors \u\ ^ than for v = 0. 



to sample the Monte Carlo averages of observables of interest. In particular 
we consider the chiral condensate, which is related to the density of the 
lowest eigenvalues of the Dirac operator by the Banks-Casher relation ()2.30|) . 
The starting point for us will be the spectral representation of the chiral 
condensate 1)2.32)1 . If we consider the dependence of the chiral condensate on 
the mass of the valence quark, m q , i.e. the quark which does not show up in 
the path integral average — then we can rewrite 

&) (m) r dz *mJ?\ Z )/V 

[mq) ~J dZ ml + {z/^V)f ' Z ~ AhV - V- 1 ) 

Now we take the first term in the expansion ()2. 100)1 and consider only the 
contribution of the lowest eigenvalue to the chiral condensate 

s&^-r* ^m, . (8.2) 



w " J« m ; + (*/(EV))= 

We evaluate this integral by a Monte Carlo method, whereby we sample 
the integral at the values z distributed according to the %RMT predictions 
for Pi'\z) given in Eq. 1)2.101)) . The resulting value can be trusted only if 
the standard deviation stabilizes. In Figure 18.61 we show the simulation of 
the integral at quark mass m q = 21.3 MeV for various topological charges. In 
particular we show the standard deviation vs. the statistics. We find that one 
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needs 0(1O 4 ) configurations to sample whereas this number reduces to 
O(10 2 ) for topological charges \u\ = 1 and 2. Therefore it is advantageous to 
simulate in topological sectors \u\ > 0. 

As a possible remedy for t he spikes there was proposed the low mode 
averaging method in Ref. |l38j |. 

8.3 Zero mode contributions to the pseudo- 
scalar correlator 

In this Section we discuss the results for the zero mode contribution to the 
pseudo-scalar correlator. The necessary theoretical introduction, as well as 
the notation, were fixed in Subsection 2.3.4. The method employed here aims 
to consider only the zero mode contributions to the pseudo-scalar correlation 
function. The latter allows to compute the low energy constants as we saw 
in Subsection 2.3.4. Here we discuss the numerical evaluation. 

We simulated quenched QCD on a 12 3 x 24 at f3 = 5.85 and on a 16 3 x 32 
at j3 — 6 lattices with periodic boundary conditions with the Wilson gauge 
action and computed the lowest eigenmodes of the Neuberger operator. On 
the 12 3 x 24 lattice at (5 = 5.85 we also calculated the lowest eigenmodes 
of the overlap hypercube Dirac operator. Out of the eigenvectors we eval- 
uated the zero mode contribution to the pseudo-scalar correlator according 
to formulae (|2.65l - l2~B7)l . We use the next-to-next-to-leading order expres- 
sion (|2.77J) in the e-regime to confront the time derivative of the zero mode 
contribution to the pseudo-scalar correlator with the corresponding quantity 
built up from the simulation data. This concept follows Ref. , which also 
presents the first numerical study. 

We looked at topological sectors \u\ — 1, 2 and 3 and we found that only 
the data from the topological sectors 1 and 2 can be used to obtain the low 
energy constants by fitting to the predictions. Figures 18771 and 18791 show 
the simulation data for the connected contribution C\ u \/L 2 and for the dis- 
connected contribution C\ v \/L 2 computed with the Neuberger operator. For 
the disconnected part we show data only for topological sectors \v\ = 2 and 3 
since the data for the topological sector 1 are the same as for the connected 
contribution up to the negative sign. In Figures 18751 we show the contribution 
of the zero modes calculated with the overlap hypercube fermions. 

The time derivative of the zero mode contribution to the pseudo-scalar 
correlator was approximated by the symmetric lattice derivative 

df /(t + 1) - /(t - 1) 

dt~ 2 • { *- 6) 
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Figure 8.7: The plot on top shows the connected contributions of zero modes 
to the pseudo-scalar correlator for topological sectors \u\ — 1, 2 and 3. The 
plot at the bottom represents the disconnected contributions to the pseudo- 
scalar correlator for topological sectors \v\ — 2 and 3. We use the Neuberger 
operator on a lattice volume V = 12 3 x 24, at f3 — 5.85. 
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Figure 8.8: The plot on top shows the connected contributions of zero modes 
to the pseudo-scalar correlator for topological sectors \v\ = 1, 2 and 3. The 
plot at the bottom represents the disconnected contributions to the pseudo- 
scalar correlator for topological sectors \v\ — 2 and 3. We use the overlap 
hypercube operator on a lattice volume V = 12 3 x 24, at /3 = 5.85. 
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Figure 8.9: The plot on top shows the connected contributions of zero modes 
to the pseudo-scalar correlator for topological sectors \v\ = 1, 2 and 3. The 
plot at the bottom represents the disconnected contributions to the pseudo- 
scalar correlator for topological sectors \v\ — 2 and 3. We use the Neuberger 
operator on a lattice volume V = 16 3 x 32, at f3 — 6. 
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First we fitted the data using the general formula (|2.66jl . Since that for- 
mula contains the parameters F n , a and (u 2 ) as unknown coefficients, we have 
to perform a 3-parameter fit. The shape coefficient (3\ for our anisotropic lat- 
tice volume was calculated according to the prescription in Ref. (ljj ]. 

Px = 0.131465 . 

The fit failed to produce reasonable estimates for the parameters. Therefore 
we followed the idea outlined in Ref. and considered the expansion of 
the formula (|2.66J) in powers of t/T — 0.5, 



D\ v \ ■ (t/T-0.5) + O((t/T-0.5) 3 
D W \ ■ (t/T- 0.5) + O((t/T-0.5f 



U) 



where the coefficients D\ v \ and D\ v \ are given by 

2\v\ f, , a fix 



Here 



D 



M 



D 



M 



M + 



2N C F2W 



+ 



I + 2 ^-2(*/ 2 >V 2 +7j7i 



2\u\ 



(F^L) 2 
13 



a 



0i 



T 2 



2N C F 2 W 



F 2 V 
+ 



r 2 \ 



5.5) 



/4(r) 

C2 



C 2 (r-0.5) + O((r 
7i(r-0.5) + 0((t 

and 7i = 

24 [ 12 



0.5) 3 ) 
0.5) 3 ) 



IE 



2^sinh^(T|pl/2) 



-0.083291(8.6) 



with |p| = (27r/L)|n|, as in Eq. ()2.78|) . The coefficients D\ v \ and D\ v \ are 
arranged as sums of terms of the leading order 0(e 4 ), next-to-leading order 
C(e 6 ) and next-to-next-to-leading order C(e 8 ), where e = l/(LF n ) and = 
F-k/V3 was introduced in Eq. (|2.7()jl . 

We employ these analytical predictions up to next-to-next-to-leading or- 
der to analyze our data. Applying one parameter fits of the data in topological 
sectors \v\ = 1, 2 and 3 to the linear part of the formula ([8.4)1 . we extracted 
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Figure 8.10: The slopes D x and D 2 vs. the fitting range [T/2 - s,T/2 + s]. 
We used the Neuberger operator on a 12 3 x 24 lattice, at /3 = 5.85 and on a 
16 3 x 32 lattice, at f3 = 6 as well as the overlap hypercube Dirac operator on 
a 12 3 x 24 lattice, at (3 = 5.85 
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the slopes of C! u <(t)/L from the data with jack-knife errors |139l Il40j . The 
results for different fitting ranges are given in Tables 18. 2| 18. 3\ 18.41 and in 
Figure EH1 

First we see that the data for the topological charge \u\ = 3 are statisti- 
cally not very stable and thus we do not use it to extract physical information. 
Note that in the case of the overlap hypercube Dirac operator the data for 
the topological sector \v\ — 3 tends to lie on the top of the data for \v\ = 2. 
Therefore we concentrate on the slopes for topological charges \v\ — 1 and 2. 

Secondly, we observe in Figure I8.1UI that the data stabilize on a fitting 
range [T/2 — s, T/2 + s] starting from s > 3. Note that we used all data points 
within this interval. Since the term of 0((t — 0.5) 3 ) in the expansion (|8.4jl is 
a NNLO term in the counting, it will be much more complicated to extract 
it. Therefore we consider only the linear term. Note that we neglect the 
0((t — 0.5) 3 ) term although we still allow the coefficients D\ v \ to be corrected 
by the NNLO terms in Eq. (fP|) . 

We now take the values for the topological susceptibility (resp. the ex- 
pectation value of the topological charge squared, (v 2 )) from our simulations 
on the 12 3 x 24 lattice at (3 — 5.85 as an input for our method, 

(y 2 ) = 10.4 ± 0.7 , for the overlap Wilson Dirac , 

(u 2 ) = 9.8 ± 0.7 , for the overlap hypercube Dirac . (8.7) 

Then we end up with a system of two equations 

D 1 = D 1 {F^a,{u 2 )) , 

D 2 = D 2 (F n ,a,(v 2 )) . (8.8) 

Here the left hand side is the numbers extracted from Figure 18.101 at s = 5 
given in Tables 18.21 - 18.41 We are aware that this choice is quite arbitrary 
but we have to postpone a more systematic analysis until larger statistics 
is available. Solving them we determine a and F w . Unlike Ref. |3(j, where 
a three parameter fit to the linear part of the expressions ()8.4j) was used to 
extract the low energy constants, in this method we employ only the values 
of the slopes and their jack-knife errors to obtain the estimates on F n and a. 
The corresponding error bars will follow as a result of the error propagation 
on the fitted values of D\ v \ and the error on the topological susceptibility, i.e. 



5a 



2 'dF v {D Xi D 2i (v*)) ;n V fdF n (D u D 2 ,(^)) 
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F n (Di, D2, (v 2 )) and ot(D\, D2, (v 2 )) are implicitly determined by the system 
of Eqs. f)8.8j) . We would like to note that the jack-knife error estimation would 
be quite difficult to apply for this method since the solutions to Eqs. ()8.8|) 
are not always real. Therefore one has to consider the resulting error bars 
with some caution. 

The results of the solution of the Eqs. ()8.8j) are given in Table 18.51 The 
first line represents F„ and a extracted from our data with the Neuberger 
operator, the values for the second line are obtained using the overlap hyper- 
cube Dirac operator. The third line represents a three parameter fit to the 
data from Ref. [3(j] (Table 18. 5|) for the Neuberger operator on lattices which 
are there denoted as C and C\. 

The results for the F v and a are compatible for different overlap fermions 
within the error bars. At this stage we only conclude that the computation 
of a poses a challenge and further investigations are needed. 

For completeness we show the statistics of our data together with the 
statistics for the C lattice in Table 18.61 

Now we fix the value of (z/ 2 ) = 10.4 for the calculations with the overlap 
Wilson Dirac operator and (u 2 ) = 9.8 for the overlap hypercube Dirac op- 
erator. Then we perform a two parameter fit of the expression (|8.4jl to the 
combined data for two topological sectors \v\ — 1 and 2 on different intervals 
t G [T/2 — s, T/2 + s\. The results for F n and a vs. the fitting range are col- 
lected in Figure IHHU We see that the values of F^ and a stabilize at s = 3. 
The values of F n and a for the Neuberger operator are consistently larger at 
finer lattice spacing. F n is still subject to renormalization which can increase 
it substantially. The value of a seems to scale within 2a. We also see that the 
results for the Neuberger operator are compatible with those of the overlap 
hypercube operator. Note, however, that the data for the overlap hypercube 
operator tend to be less noisy than for the Neuberger operator on the lattice 
12 3 x 24, at = 5.85. 

The estimates and the jack-knife error bars of F n and a, as well for the 
re-calculated D\ u \, for the fitting range [T/2 — 5, T/2 + 5] are collected in 
Table 18.71 First of all we can confirm again that the given statistics is not 
sufficient to obtain a stable value of the slope D% at least on 12 3 x 24 lattice 
at (3 — 5.85. Next we see that the two methods that we discussed lead 
to compatible results. However the two parameter fit produces significantly 
more reliable error bars. Note that if we change the topological susceptibility 
in our two parameter fit by 10% this will lead only to a 0(10%) change in 
a while F^ remains almost the same. We conclude that performing the two 
parameter fit is still the best method to obtain the estimates for F„ and a with 
a control over the errors. We argue that the calculation of the mere zero mode 
contribution is potentially capable of the extraction of F n , a and topological 
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Figure 8.11: F n and a extracted from a two parameter fit to the combined 
data from topological sectors \v\ = 1 and 2 vs. the fitting range [T/2 — 
s,T/2 + s]. We used the Neuberger operator on a 12 3 x 24 lattice at /3 — 5.85 
and on a 16 3 x 32 lattice at (3 = 6, as well as the overlap hypercube Dirac 
operator on a 12 3 x 24 lattice, at (3 = 5.85. 
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susceptibility Xt once a good statistics is available. Its obvious virtue is that it 
is not affected by the problems of small quark masses and small eigenvalues of 
the Dirac operator which magnify the noise of the propagator as we discussed 
in Section 8.2. 
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lattice 
size 




n 


physical 
volume 


numbe 
i/ = 


r of confi 
\u\ = 1 


gurations 
\v\ = 2 


12 4 


6 


1.4 


(1.12 fm) 4 


47 


78 


24 


10 3 x 24 


6 


1.4 


(0.93 fm) 3 x 1.86 fm 


20 


24 


17 



Table 8.1: The statistics of our simulations on two lattice sizes. 



fitting 
range 


H = i 


H = 2 


D\»\ 
\v\ = 3 


X 2 /d.o.f. 


X 2 /d.o.f. 


X 2 /d.o.f. 


[11,13] 


9.2 ±4.5 


22.7 ±8.6 


20.3 ± 12.7 


0.2 


0.05 


4.9 


[10, 14] 


9.7 ±3.8 


19.7 ±6.2 


24.3 ±9.9 


0.1 


0.05 


2.8 


[7, 17] 


12.2 ± 1.6 


20.2 ±3.0 


27.0 ±±5.0 


1.3 


0.1 


1.3 



Table 8.2: The slopes D\ v \ vs. the fitting range computed with the overlap 
Wilson Dirac operator at (3 — 5.85 on 12 3 x 24 lattice. 



fitting 
range 


\v\ = 1 


\u\ = 2 


\v\ = 3 


X 2 /d.o.f. 


X 2 /d.o.f. 


X 2 /d.o.f. 


[11,13] 


19.4 ±3.3 


30.6 ±5.7 


27.7 ± 11.4 


0.16 


1.6 


0.05 


[10, 14] 


17.6 ±2.6 


27.1 ±4.4 


23.7 ±8.9 


0.17 


0.86 


0.07 


[7, 17] 


14.2 ±1.2 


26.2 ±2.7 


26.2 ±4.5 


0.6 


0.4 


0.15 



Table 8.3: The slopes D\ v \ vs. the fitting range computed with the overlap 
hypercube Dirac operator at (3 = 5.85 on 12 3 x 24 lattice. 



fitting 
range 


\v\ = 1 


\u\ = 2 


D\„ 

\u\ = 3 


X 2 /d.o.f. 


X 2 /d.o.f. 


X 2 /d.o.f. 


[15, 17] 


3.5 ±7.8 


5.1 ± 11.8 


48.1 ± 14.5 


0.21 


0.87 


0.47 


[14, 18] 


4.9 ±5.4 


14.7 ±8.8 


35.4 ± 12.9 


0.12 


0.58 


0.32 


[11,21] 


8.2 ±2.7 


25.0 ±5.5 


22.0 ±7.4 


0.44 


0.47 


0.3 



Table 8.4: The slopes D\ v \ vs. the fitting range computed with the overlap 
Wilson Dirac operator at (3 = 6 on 16 3 x 32 lattice. 



p 


V 


D x 


D 2 


F n [MeV] 


a 


remarks 


5.85 


12 3 x 24 


13.5 ± 1.4 


21.8 ±2.6 


80.7 ±9.1 


-12.6 ±6.6 


Neuberger 


5.85 


12 3 x 24 


14.2 ±1.2 


26.2 ±2.7 


85.7 ± 11 


-6.2 ±8.2 


overlap HF 


6 


20 4 


3.9 ±0.5 


9.4 ±0.7 






Ci 


5.8784 


16 4 


3.4 ±0.3 


8.9 ±0.5 


(80.8,107.4) 


(-1.8,7.8) 


Co 



Table 8.5: The first two lines show F n and a as a solution of the Eq. (|8.5j) . 
extracted from our data. The third line gives for comparison the values quoted 
in Ref. 0. 
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operator 
kernel 


lattice 
size 


P 


physical 
volume 


number 

H = i 


of config 
H = 2 


urations 
\v\ = 3 


Wilson 


12 3 x 24 


5.85 


(1.48 fm) 3 x 2.96 fm 


90 


72 


64 


Wilson 


16 3 x 32 


6 


(1.49 fm) 3 x 2.98 fm 


92 


70 


71 


hypercube 


12 3 x 24 


5.85 


(1.48 fm) 3 x 2.96 fm 


91 


83 


69 


Wilson 


16 4 


5.8784 


(1.86 fm) 4 


229 


186 





Table 8.6: The statistics for our data (1st, 2nd and 3rd lines) and for the Co 
lattice in Ref. (fourth line). 



Di 


D 2 


D 3 




a 


rVd.o.f. 


12.4 ± 1.6 


20.2 ±3.1 


7 A ±17 


79.9±8.7MeV 


-13.0 ±6.1 


0.7 


14.3 ± 1.3 


26.2 ±2.7 


23.3 ± 15 


85.1 ± 10.5 MeV 


-3.9 ±7.0 


0.5 


9.9 ± 1.9 


21.8 ±3.8 


33.7 ±5.7 


132.9 ± 0.0 MeV 


16.5 ±5.7 


0.6 



Table 8.7: Results of the two parameter fit in the range [T/2 — 5,T/2 ± 5] 
for the combined data in topological sectors \v\ — 1 and 2. The first line 
represents the results for the overlap Wilson Dirac operator and the second 
line represents the results for overlap hypercube Dirac operator on a lattice 
of size 12 3 x 24, at /3 = 5. 85. The third line represents the results for the 
overlap Wilson Dirac operator on a lattice of size 16 3 x 32, at (3 = 6. 
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Topology conserving Liischer 
gauge action 

9.1 Motivation 

We employ the overlap operator as a regularization for the fermion fields. 
It obeys the Ginsparg- Wilson relation (|2.122jl and thus the chiral symme- 
try is preserved at a finite lattice spacing as discussed in Ref. [69] and we 
demonstrated in Section 2.5.3. This enables us to perform simulations at 
small quark masses and to define the topological charge as the index of the 
overlap operator. In particular it is applied to explore the e-regime. In the 
previous Chapter we saw that in order to perform the precision calculation of 
the low energy constants in the e-regime one has to collect as much statistics 
as possible for a given topological sector. A computer time consuming part 
of these calculations is that for each configuration one needs to compute the 
index of the overlap operator. In particular we saw that the simulations in 
the topological sector v — 0, which occurs quite often in the simulations, 
are plagued by the large statistical fluctuations and therefore it is favorable 
to evaluate observables in the sectors of \v\ ^ 0. As we learned from the 
xRMT study of the probability distributions of the individual eigenvalues, 
the topological sectors \u\ — 1 and 2 are of primary interest for the simu- 
lations in moderate physical volumes. In simulations on a large lattices the 
charges fluctuate between v — — 10 and 10, see Figure IT7T1 Therefore it would 
be highly profitable for the simulations if a method was found to sample the 
statistics within a certain topological sector without recomputing the index. 
In other words, one seeks a method which provides a large autocorrelation 
time for the topological charge and still small autocorrelation times for the 
observables of interest. 
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The change of the topological sector for the overlap Dirac operator implies 
that an eigenvalue flips from the origin of the Ginsparg- Wilson circle to 2/x/a, 
i.e. to the cut-off scale, or vice versa as the gauge field is deformed. If the 
corresponding deformation is done continuously then the inverse square root 
of A' A will develop a singularity. As it was already mentioned it is guaranteed 
under the con ditio n ()3.19j) that A* A remains positive. It was shown later by 
H. Neuberger |l4l| that this condition can be relaxed to 

e < ~ _L (Q l) 

~ + 20.5 • V ' } 

If it is fulfilled then the topological charge is guaranteed not to change under 
a continuous deformation of the gauge field. A non-continuous deformation 
can, however, still change it but the hope is that the probability of such a 
change is small in simulations. Applying the weak coupling expansion for the 
Wilson gauge action we found that the above relation is satisfied for the mean 
plaquette value at (3 ~ 40. Obviously simulations are not possible with such 
a large value of (3 due to the strong autocorrelations and tiny lattice spacings. 
However one hopes that the topological charge will be still conserved over 
a large number of configurations even if the gauge action only favors the 
plaquettes obeying the further relaxed condition (JHHJ). 

It is possible to introduce a lattice gauge action which will suppress large 
values of Sp(Up) = 1 — (l/3)ReTr£/p and thereby the condition f)3. 18j) will 
be fulfilled. This will in turn lead to conservation of the topological charge 
along the simulation trajectory. To this end we define 

I Sp{uj P Jy) q (tt{P) \ ^ 

oo S P (UiJu) > e 



where a > 0. In particular, a = 1 was used by M. Lus cher |142l |143] for 



conceptual purposes and was also applied in Ref. 1144| to the Schwinger 
model. In the latter the authors used e = 1 as opposed to the strict theoretical 
bound of 1/5 for two dimensions, in order to let the gauge fields fluctuate 
more while still suppressing topological transitions. This was sufficient to 
stabilize the topological history over hundreds of configurations. A potential 
problem for the simulations with this action can be that due to the stringent 
constraint the resulting lattice spacing may become too small. This would 
require very large lattice volumes to obtain a reasonable physical volume for 
simulations. Yet another caveat was raised by M. Creutz that such an action 
may lead to a non-positive transfer matrix which wil l ma ke the transition to 
Minkowski space at finite lattice spacing impossible 145l |. However, one can 



think of first taking the continuum limit and then going to Minkowski space. 
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9.2 Results on the topological histories 

We used the local Hybrid Monte Carlo (HMC) algorithm — for a descrip- 
tion we refer to Appendix A — to simulate the action (|9.2j) with a = 1 in 
quenched QCD. The standard tools like heat bath and overrelaxation are 
not applicable in this case since they very much rely on the properties of the 
gauge action (|H.2jl which do not hold for the action (|9.2jl . For HMC one has 
to compute the "force" F x>a , 

sSqju^L) _ Wilson i + ^-Sp 



5U X ,„ " x ^ a (l-S P /e)^ 



(P) 

pWilson vl->P\vx,hvj ,„ „s 

Figure iniU shows the dependence of S a /Sp vs. Sp as well as F a /F P vs. 
Sp. As a decreases the "force" becomes more flat in the domain of Sp close 
to the free field. 

If the new action is defined with quite a stringent constraint this will 
reduce the freedom for the gauge field to fluctuate which could imply a very 
small lattice spacing as we mentioned before. To investigate this we start 
off by searching a set of values of {e, (3} which will provide the same (Sp) 
in simulations on 4 4 lattice. The obtained curve of the "constant physics" 
was used as a guideline for the further analysis where we computed r$ for 
the specific values of e, j3 on the curve. The results for the curve of the 
"constant physics" are reported in Figure 19.21 As e decreases the value of j3 
must be also decreased in order to obtain the same expectation value for the 
plaquette (Sp). Another problem is the rejection of a link variable due to the 
violation of the boundary in definition of the action (|9.2jl . It t urned out that 
the rejection rate is increased as (3 is decreased, see Table EH Ei. 

Next we picked out a few points from the curve and computed r on a 
16 4 lattice. The results for r are given in Table l9~Tl We see that r slightly 
increases as j3 is decreased. As a pilot test for the topological stability we em- 
ployed the cooling method |l47| . This is a relatively fast way compared to the 
overlap index to estimate the topological charge of a gauge configuration at 
moderate (3. In this approach one deforms a gauge configuration by minimiz- 
ing its gauge action locally to reach the approximate solutions of the lattice 
Yang-Mills equations of motion, which are the Euler equations for the Wil- 
son Yang-Mills gauge action. The resulting gauge configurations, which are 
picked out from the first stable plateau, dissociate in objects with topological 
charges v — ±1. They can be interpreted as instantons and anti-instantons. 
Since the cooled configurations are smooth, these charges can be determined 



5Sp(U : 
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Figure 9.1: The plot on top shows the dependence of S a /Sp vs. Sp. The 
plot at the bottom shows the dependence of the "force" F a / F P vs. Sp. The 
theoretical bound is shown as an asymptotic vertical line at Sp — 1/e. 
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Figure 9.2: The lines of constant expectation value of the plaquette, (Sp), in 
the plane spanned by the action parameters 1/e and (3. 



using a simple la ttice discretization of the continuum formula (|2.8|) given for 
instance in Ref. |l48j |. Hence the topological charge of a cooled gauge con- 
figuration can be determined by the difference of number of instantons iV + 
and anti- instantons iV_. To quantify the topological charge conservation we 
computed the autocorrelation time of Q C ooi = N + — iV_. The main conjecture 
is that the autocorrelation time of this quantity is a good estimate for the 
autocorrelation time of the topological charge given by the index of the over- 
lap operator for the equilibrium configurations. The latter is supported by 
the fact that the corresponding equilibrium configurations are quite smooth 
(@w ~ 6.19). Hence the definition of the index given by cooling match es th e 
definition of the index given by the overlap operator at least up to 80% 149l | . 



This is verified by our study at stronger coupling (3 = 6, shown in Figure! 
It shows the distribution of the difference of the topological charge defined 
with cooling and the topological charge defined by the Neuberger index. 

In Figures 19.41 we show the histories of Q C ooi for the parameters (3 and e 
given in Table I9~T1 

The autocorrelation times for the quantity Q coo \, t^ 1 , as well as for the 
plaquette Sp, r^ q , are given in Table I9~T1 We see that r^ 1 becomes larger 
as e is decreased, while the autocorrelation times for the plaquette become 
even smaller. This gives confidence that the pair 1/e = 1.52 and (3 = 0.3 
constitutes a good parameter set for the action to conduct simulations in the 
e-regime. 
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1 1 1 1 1 1 1 1 1 

-2-10 1 2 

Figure 9.3: The histogram of the difference between topological charge defined 
by the cooling, Q coo \, and the topological charge defined by the index of the 
Neuberger operator, (^Neuberger- For our parameters V = 10 3 x 24, (3 = 6 we 
observe 80 % agreement. 





P 


rejection rate 


r /a 


Pw 


_cool 
aut 


plaq 
' aut 





6.19 





7.25 


6.19 


3.94 


7.27 


1.25 


0.8 


< 0.1% 


7.0(1) 


6.17 


5.07 


1.15 


1.52 


0.3 


5% 


7.3(4) 


6.19 


21.03 


0.85 



Table 9.1: Results for the parameter r /a computed in simulations on 16 4 
lattice for various values of e and (3. We give the corresponding value of (3 
for the Wilson action, which we denote as /3w The autocorrelation times for 
Qcooi an d for the plaquette are reported. We also show the the rejection rate 
for the links. 
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Figure 9.4: Histories of Q cool at {/3,1/e} = {6.19,0}, {0.8,1.25}, {0.3,1.52} 
for the plot on top, in the middle and at the bottom, respectively. 
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In this thesis we presented a study of the feasibility to extract physical infor- 
mation in the e-regime using various formulations of chiral lattice fermions. 
In this final Section we briefly summarize the main results. 

The investigations were done in the framework of quenched lattice QCD. 
To discretize the fermion fields we chose the Ginsparg- Wilson fermions. In 
particular we applied the Neuberger operator and the overlap hypercube op- 
erator. These formulations provide a regularization scheme which preserves 
the chiral symmetry at a finite lattice spacing. As a consequence the topolog- 
ical charge is well defined by the index of the overlap Dirac operator. There- 
fore by simulating quenched QCD in the e-regime one can reach the domain 
of very small quark masses. The computations are, however, quite tedious 
and therefore at present only the quenched approximation can be studied. 
Controlling the finite size effects in QCD simulations in the e-regime one 
can obtain the low energy constants of the chiral Lagrangian, which are also 
those of the infinite volume up to quenching effects. 

One of the main topics of these investigations was the analysis of the 
probability distributions of the individual eigenvalues of the Dirac operator 
in the QCD simulations in the e-regime. We confronted the simulation data 
with the analytical prediction conjectured in the framework of the xRMT. 
The conjecture is that the xRMT is a theory which can describe eigenvalue 
densities in the e-regime. It is mainly based on the fact that in the limit of 
large matrices the xRMT partition function coincides with the xPT parti- 
tion function to the leading order in the e-regime. Also some quantities like 
the microscopic spectral density were calculated for the special cases both, in 
xRMT and in xPT, leading to the same expressions. The ability of xR-MT 
to explain the simulations in the e-regime being yet a conjecture received 
a convincing support in our investigations. To this end we found that the 
behavior of the probability distributions of the individual eigenvalues can be 
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reliably described by the xRMT if the box length of the physical volume ex- 
ceeds about L > 1.12 fm. The chiral condensate enters as the only parameter 
for the probability distributions of the individual eigenvalues. Therefore it 
can be obtained by fitting the xRMT curves to the simulation data. 

A next result obtained from this study was the determination of the lower 
bound for the size of the physical volume where the data can be described 
by %PT. It has proven to help in understanding the behavior of the two- 
point mesonic functions. In particular we could explain why the correlators 
calculated in quenched simulations on (0.93 fm) 3 x 1.86 fm are not modeled 
by the x?T. It turned out that the size of the physical volume of L as 0.93 fm 
is too small for the simulations in the e-regime and a size of L m 1.12 fm 
lies on the threshold. When trying to extract the low energy constants F n 
and E from the axial- vector two-point function we found that F n , which 
enters the formula in an additive constant, can be obtained by fitting. The 
chiral condensate E, on the contrary, enters in a way that leaves the analytic 
formula almost insensitive to drastic changes in E for the topological sectors 
\u\ > 0. The sensitivity is restored only for v — 0. This makes it prohibitively 
difficult to extract E using the axial- vector two-point function at \v\ ^ 0. On 
the other hand, the signal for the axial-vector correlator in the topological 
sector v — turned out to be very noisy. This we attributed to a high 
probability of occurrence of near- zero modes. The latter was also supported 
by our xRMT study of the eigenvalues. 

To estimate the statistics needed for the topological sectors \v\ — 0, 1 
and 2 we considered the contribution to the chiral condensate given by the 
distribution of the lowest eigenvalue of the Dirac operator E^ n . Then we sim- 
ulated fake eigenvalues with these distributions and obtained the necessary 
statistics to stabilize the standard deviation of the E^ n . It turned out that 
one needs a tremendous statistics of O(10 4 ) configurations in the topological 
sector v — (for this stabilization). For higher topological sectors i\v\ = 1, 
2) this number is decreased to O(10 2 ) configurations. 

As an alternative method to calculate the low energy constants one can 
consider only the contribution of the zero modes to the mesonic two-point 
functions. We computed such contributions to the pseudo-scalar two-point 
function. The analytical expressions for them are known and we fitted them 
to our data to explore the possibility of extracting the parameters F n and 
a — yet another low energy parameter appearing in the chiral quenched 
Lagrangian. We were able to obtain reasonable bounds on these values with 
our statistics. Therefore we conclude that this is a promising method to 
compute the low energy constants. 

Next we considered the topological susceptibility. This quantity was com- 
puted using two formulations of the overlap operator: the Neuberger fermions 
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(with the usual Wilson action as its kernel) and the hypercube fermions (with 
the truncated fixed point action as its kernel). The overlap hypercube opera- 
tor is expected to be superior to the conventional formulation. It turned out 
that the locality is improved and the topological susceptibility comes closer 
to the continuum value at the same lattice spacing. 

As a theoretical development of the simulations in the e-regime one can 
think of a method to simulate a set of configurations belonging to a certain 
topological charge without recomputing the index. The problem can be re- 
formulated as to simulate configurations with the topological charge being 
strongly autocorrelated while other observables of interest remaining decor- 
related. This we studied using a gauge action proposed by M. Liischer. We 
are confident that we found a good set of parameters for the Liischer gauge 
action to meet these requirements. 
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Appendix A 

A local Hybrid Monte Carlo 
algorithm 



In this Ap pendix we discuss the implementation of the local HMC algo- 
rithm 150l ] that we used to simulate th e Luscher gauge action. It is based on 



the HMC algorithm described in Refs. [151, 15^ 

Let us consider the Hamiltonian of the 577(3) gauge theory 



1 8 



r- 



We define ir x ^ = Y^j=i ^x^j^j, where are the Gell-Mann matrices for the 
SU(3) gauge group. The variables n x ^ are canonical momenta for U Xsft . The 
local HMC algorithm for the simulation of the gauge fields U Xtfl G SU(3) 
starts off by choosing a random value of the momentum n X)fl according to 
the Gaussian distribution 

pgka « ri ex p {-^'.wj • ( a - 2 ) 

j=l < > 

Then, using the Hamiltonian H the pair [ir x ,n, U x ^} is changed along 
the discretized trajectory in phase space up to the end point [7r' , U' ] = 
T H [ir,U], where T H is the discretized evolution operator corresponding to 
the Hamiltonian H. It is assumed that the discretized classical dynam- 
ics equations, which define the operator Tjj, are such that the mapping 
[ n x,tM, U x>l t] — > [ir', U' ] is reversible in the sense 

Ph(K», E£J - Km, U X J) = Ph{[-^ U x J <- U' X J) . (A.3) 
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Ph is the probability of the configuration transition 



Ph(K,^ U' x J «- [tt^, ^J) = <5([tt^, t^J - T^, C^J) , (A.4) 

which takes this form since the change is deterministic in the case of the 
classical equations of motion. The reversibility condition holds for instance 
for the leap frog integration. After the trajectory [ir Xtfl , U XtfM ] — > [n^u, U'] is 
determined, the pair [k x ^, U' x ^} is either accepted or rejected, and the old pair 
[^x,/i7 U Xtfl ] is kept. The acceptance is defined by the Metropolis probability 

Pa(K^U' x J - [vr^,^) =min{l,e- A ^[<^J+AHK M ,^ M ]| j (A5) 
where 

8 1 

Ai? [tt w u x j = J2 + E p spWSL) ( a - 6 ) 

i=i 

is the part of the Hamiltonian H which includes tt x ^ and U x ^. The total 
transition probability is 

Ph(K^ U' x J «- K„ 17.J) P g [tt x J . (A.7) 

It can be shown that this probability distribution has the canonical distribu- 
tion proportional to exp (—Sym[U]) as a unique fixed point. It can be proved 
by using the reversibility condition (|A.3|) . Therefore a finite difference ap- 
proximation to the Hamiltonian equations must also account for it. 

Applying now this algorithm to each pair [vr^^, U XjfM ] we perform the up- 
date of the gauge field throughout the lattice. 

We close the explanation of the algorithm by giving the description of 
the finite difference approximation to the Hamiltonian equations, i.e. to the 
evolution operator Th- It is given by the leapfrog integration, which proceeds 
as follows: one performs a first half step 

ir x ,»,j{— ) = n Xilhj {0) jjj ' ^ 8 ^ 

3 x ^ 

and then for steps k — 0, 1, . . . , n — 1 
U X) ^{kAr + At) = exp j^zA^^Ar + 4r) Ar j U x ^(kAr) , (A.9) 



and for k — 1, 2, . . . , n — 1 

^(JfcAr + = W *Ar - - ^^ A T , (A.10) 
and finally 

= - - iA ^^'" AT)l Ar . (A.11) 

where ASym[U x ,h] = X^<^ P Sp(U x P J v ) is the part of the lattice gauge action 
which contains U XtjX . The derivative appearing in these formulae is defined as 

W^jj - A f \ e **iTj 1 f A 12) 
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